Changes in trunk/python/scantable.py [2687:3010]

File:

• trunk/python/scantable.py (modified) (81 diffs)

trunk/python/scantable.py

@@ -2 +2 @@
 
 import os
+import re
 import tempfile
 import numpy
@@ -45 +46 @@
         l=f.readline()
         f.close()
-        if ( l.find('Scantable') != -1 ):
-            return True
-        elif ( l.find('Measurement Set') == -1 and 
-               l.find('Image') == -1 ):
+        match_pattern = '^Type = (Scantable)? *$'
+        if re.match(match_pattern,l):
             return True
         else:
@@ -175 +174 @@
     return str(showprogress).lower() + ',' + str(minnrow)
 
+def pack_blinfo(blinfo, maxirow):
+    """\
+    convert a dictionary or a list of dictionaries of baseline info
+    into a list of comma-separated strings. 
+    """
+    if isinstance(blinfo, dict):
+        res = do_pack_blinfo(blinfo, maxirow)
+        return [res] if res != '' else []
+    elif isinstance(blinfo, list) or isinstance(blinfo, tuple):
+        res = []
+        for i in xrange(len(blinfo)):
+            resi = do_pack_blinfo(blinfo[i], maxirow)
+            if resi != '':
+                res.append(resi)
+    return res
+
+def do_pack_blinfo(blinfo, maxirow):
+    """\
+    convert a dictionary of baseline info for a spectrum into
+    a comma-separated string. 
+    """
+    dinfo = {}
+    for key in ['row', 'blfunc', 'masklist']:
+        if blinfo.has_key(key):
+            val = blinfo[key]
+            if key == 'row':
+                irow = val
+            if isinstance(val, list) or isinstance(val, tuple):
+                slval = []
+                for i in xrange(len(val)):
+                    if isinstance(val[i], list) or isinstance(val[i], tuple):
+                        for j in xrange(len(val[i])):
+                            slval.append(str(val[i][j]))
+                    else:
+                        slval.append(str(val[i]))
+                sval = ",".join(slval)
+            else:
+                sval = str(val)
+
+            dinfo[key] = sval
+        else:
+            raise ValueError("'"+key+"' is missing in blinfo.")
+
+    if irow >= maxirow: return ''
+    
+    for key in ['order', 'npiece', 'nwave']:
+        if blinfo.has_key(key):
+            val = blinfo[key]
+            if isinstance(val, list) or isinstance(val, tuple):
+                slval = []
+                for i in xrange(len(val)):
+                    slval.append(str(val[i]))
+                sval = ",".join(slval)
+            else:
+                sval = str(val)
+            dinfo[key] = sval
+
+    fspec_keys = {'poly': 'order', 'chebyshev': 'order', 'cspline': 'npiece', 'sinusoid': 'nwave'}
+
+    fspec_key = fspec_keys[dinfo['blfunc']]
+    if not blinfo.has_key(fspec_key):
+        raise ValueError("'"+fspec_key+"' is missing in blinfo.")
+
+    clip_params_n = 0
+    for key in ['clipthresh', 'clipniter']:
+        if blinfo.has_key(key):
+            clip_params_n += 1
+            dinfo[key] = str(blinfo[key])
+
+    if clip_params_n == 0:
+        dinfo['clipthresh'] = '0.0'
+        dinfo['clipniter']  = '0'
+    elif clip_params_n != 2:
+        raise ValueError("both 'clipthresh' and 'clipniter' must be given for n-sigma clipping.")
+
+    lf_params_n = 0
+    for key in ['thresh', 'edge', 'chan_avg_limit']:
+        if blinfo.has_key(key):
+            lf_params_n += 1
+            val = blinfo[key]
+            if isinstance(val, list) or isinstance(val, tuple):
+                slval = []
+                for i in xrange(len(val)):
+                    slval.append(str(val[i]))
+                sval = ",".join(slval)
+            else:
+                sval = str(val)
+            dinfo[key] = sval
+
+    if lf_params_n == 3:
+        dinfo['use_linefinder'] = 'true'
+    elif lf_params_n == 0:
+        dinfo['use_linefinder'] = 'false'
+        dinfo['thresh']         = ''
+        dinfo['edge']           = ''
+        dinfo['chan_avg_limit'] = ''
+    else:
+        raise ValueError("all of 'thresh', 'edge' and 'chan_avg_limit' must be given to use linefinder.")
+    
+    slblinfo = [dinfo['row'], dinfo['blfunc'], dinfo[fspec_key], dinfo['masklist'], \
+                dinfo['clipthresh'], dinfo['clipniter'], \
+                dinfo['use_linefinder'], dinfo['thresh'], dinfo['edge'], dinfo['chan_avg_limit']]
+    
+    return ":".join(slblinfo)
+
+def parse_fitresult(sres):
+    """\
+    Parse the returned value of apply_bltable() or sub_baseline() and
+    extract row number, the best-fit coefficients and rms, then pack
+    them into a dictionary.
+    The input value is generated by Scantable::packFittingResults() and
+    formatted as 'row:coeff[0],coeff[1],..,coeff[n-1]:rms'. 
+    """
+    res = []
+    for i in xrange(len(sres)):
+        (srow, scoeff, srms) = sres[i].split(":")
+        row = int(srow)
+        rms = float(srms)
+        lscoeff = scoeff.split(",")
+        coeff = []
+        for j in xrange(len(lscoeff)):
+            coeff.append(float(lscoeff[j]))
+        res.append({'row': row, 'coeff': coeff, 'rms': rms})
+
+    return res
+
+def is_number(s):
+    s = s.strip()
+    res = True
+    try:
+        a = float(s)
+        res = True
+    except:
+        res = False
+    finally:
+        return res
+
+def is_frequency(s):
+    s = s.strip()
+    return (s[-2:].lower() == "hz")
+
+def get_freq_by_string(s1, s2):
+    if not (is_number(s1) and is_frequency(s2)):
+        raise RuntimeError("Invalid input string.")
+    
+    prefix_list = ["a", "f", "p", "n", "u", "m", ".", "k", "M", "G", "T", "P", "E"]
+    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]
+
+    s1 = s1.strip()
+    s2 = s2.strip()
+    
+    prefix = s2[-3:-2]
+    if is_number(prefix):
+        res1 = float(s1)
+        res2 = float(s2[:-2])
+    else:
+        factor = factor_list[prefix_list.index(prefix)]
+        res1 = float(s1) * factor
+        res2 = float(s2[:-3]) * factor
+
+    return (res1, res2)
+
+def is_velocity(s):
+    s = s.strip()
+    return (s[-3:].lower() == "m/s")
+
+def get_velocity_by_string(s1, s2):
+    if not (is_number(s1) and is_velocity(s2)):
+        raise RuntimeError("Invalid input string.")
+
+    # note that the default velocity unit is km/s
+    prefix_list = [".", "k"]
+    factor_list = [1e-3, 1.0]
+
+    s1 = s1.strip()
+    s2 = s2.strip()
+
+    prefix = s2[-4:-3]
+    if is_number(prefix): # in case velocity unit m/s
+        res1 = float(s1) * 1e-3
+        res2 = float(s2[:-3]) * 1e-3
+    else:
+        factor = factor_list[prefix_list.index(prefix)]
+        res1 = float(s1) * factor
+        res2 = float(s2[:-4]) * factor
+
+    return (res1, res2)
+
+def get_frequency_by_velocity(restfreq, vel, doppler):
+    # vel is in unit of km/s
+
+    # speed of light
+    vel_c = 299792.458
+
+    import math
+    r = vel / vel_c
+
+    if doppler.lower() == 'radio':
+        return restfreq * (1.0 - r)
+    if doppler.lower() == 'optical':
+        return restfreq / (1.0 + r)
+    else:
+        return restfreq * math.sqrt((1.0 - r) / (1.0 + r))
+
+def get_restfreq_in_Hz(s_restfreq):
+    value = 0.0
+    unit = ""
+    s = s_restfreq.replace(" ","")
+
+    for i in range(len(s))[::-1]:
+        if s[i].isalpha():
+            unit = s[i] + unit
+        else:
+            value = float(s[0:i+1])
+            break
+
+    if (unit == "") or (unit.lower() == "hz"):
+        return value
+    elif (len(unit) == 3) and (unit[1:3].lower() == "hz"):
+        unitprefix = unit[0]
+        factor = 1.0
+
+        prefix_list = ["a", "f", "p", "n", "u", "m", ".", "k", "M", "G", "T", "P", "E"]
+        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]
+        factor = factor_list[prefix_list.index(unitprefix)]
+        """
+        if (unitprefix == 'a'):
+            factor = 1.0e-18
+        elif (unitprefix == 'f'):
+            factor = 1.0e-15
+        elif (unitprefix == 'p'):
+            factor = 1.0e-12
+        elif (unitprefix == 'n'):
+            factor = 1.0e-9
+        elif (unitprefix == 'u'):
+            factor = 1.0e-6
+        elif (unitprefix == 'm'):
+            factor = 1.0e-3
+        elif (unitprefix == 'k'):
+            factor = 1.0e+3
+        elif (unitprefix == 'M'):
+            factor = 1.0e+6
+        elif (unitprefix == 'G'):
+            factor = 1.0e+9
+        elif (unitprefix == 'T'):
+            factor = 1.0e+12
+        elif (unitprefix == 'P'):
+            factor = 1.0e+15
+        elif (unitprefix == 'E'):
+            factor = 1.0e+18
+        """
+        return value*factor
+    else:
+        mesg = "wrong unit of restfreq."
+        raise Exception, mesg
+
+def normalise_restfreq(in_restfreq):
+    if isinstance(in_restfreq, float):
+        return in_restfreq
+    elif isinstance(in_restfreq, int) or isinstance(in_restfreq, long):
+        return float(in_restfreq)
+    elif isinstance(in_restfreq, str):
+        return get_restfreq_in_Hz(in_restfreq)
+    elif isinstance(in_restfreq, list) or isinstance(in_restfreq, numpy.ndarray):
+        if isinstance(in_restfreq, numpy.ndarray):
+            if len(in_restfreq.shape) > 1:
+                mesg = "given in numpy.ndarray, in_restfreq must be 1-D."
+                raise Exception, mesg
+        
+        res = []
+        for i in xrange(len(in_restfreq)):
+            elem = in_restfreq[i]
+            if isinstance(elem, float):
+                res.append(elem)
+            elif isinstance(elem, int) or isinstance(elem, long):
+                res.append(float(elem))
+            elif isinstance(elem, str):
+                res.append(get_restfreq_in_Hz(elem))
+            elif isinstance(elem, dict):
+                if isinstance(elem["value"], float):
+                    res.append(elem)
+                elif isinstance(elem["value"], int):
+                    dictelem = {}
+                    dictelem["name"]  = elem["name"]
+                    dictelem["value"] = float(elem["value"])
+                    res.append(dictelem)
+                elif isinstance(elem["value"], str):
+                    dictelem = {}
+                    dictelem["name"]  = elem["name"]
+                    dictelem["value"] = get_restfreq_in_Hz(elem["value"])
+                    res.append(dictelem)
+            else:
+                mesg = "restfreq elements must be float, int, or string."
+                raise Exception, mesg
+        return res
+    else:
+        mesg = "wrong type of restfreq given."
+        raise Exception, mesg
+
+def set_restfreq(s, restfreq):
+    rfset = (restfreq != '') and (restfreq != [])
+    if rfset:
+        s.set_restfreqs(normalise_restfreq(restfreq))
+
 class scantable(Scantable):
     """\
@@ -210 +513 @@
                           Default (false) is taken from rc file.
 
+            getpt:        Whether to import direction from MS/POINTING
+                          table properly or not.
+                          This is effective only when filename is MS.
+                          The default (True) is to import direction
+                          from MS/POINTING.
         """
         if average is None:
@@ -244 +552 @@
                     self._fill([filename], unit, average, opts)
                 elif os.path.isfile(filename):
-                    self._fill([filename], unit, average)
+                    opts={'nro': {}}
+                    nrokeys=['freqref']
+                    for key in nrokeys:
+                        if key in args.keys():
+                            opts['nro'][key] = args[key]
+                    self._fill([filename], unit, average, opts)
                     # only apply to new data not "copy constructor"
                     self.parallactify(parallactify)
@@ -578 +891 @@
 
     @asaplog_post_dec
-    def stats(self, stat='stddev', mask=None, form='3.3f', row=None):
+    def stats(self, stat='stddev', mask=None, form='3.3f', row=None, skip_flaggedrow=False):
         """\
         Determine the specified statistic of the current beam/if/pol
@@ -596 +909 @@
             row:     row number of spectrum to process.
                      (default is None: for all rows)
+
+            skip_flaggedrow: if True, skip outputting text for flagged
+                             spectra. default is False.
 
         Example:
@@ -608 +924 @@
                              "number of channels. Please use setselection() "
                              "to select individual IFs")
-        rtnabc = False
-        if stat.lower().endswith('_abc'): rtnabc = True
         getchan = False
         if stat.lower().startswith('min') or stat.lower().startswith('max'):
@@ -615 +929 @@
             getchan = True
             statvals = []
-        if not rtnabc:
+
+        rtnabc = False
+        if stat.lower().endswith('_abc'):
+            rtnabc = True
+        else:
             if row == None:
                 statvals = self._math._stats(self, mask, stat)
@@ -641 +959 @@
             statunit= ''
             if getchan:
-                qx, qy = self.chan2data(rowno=i, chan=chan[i])
-                if rtnabc:
-                    statvals.append(qx['value'])
-                    refstr = ('(value: %'+form) % (qy['value'])+' ['+qy['unit']+'])'
-                    statunit= '['+qx['unit']+']'
+                if self._is_all_chan_flagged(i):
+                    if rtnabc:
+                        statvals.append(None)
                 else:
-                    refstr = ('(@ %'+form) % (qx['value'])+' ['+qx['unit']+'])'
+                    qx, qy = self.chan2data(rowno=i, chan=chan[i])
+                    if rtnabc:
+                        statvals.append(qx['value'])
+                        refstr = ('(value: %'+form) % (qy['value'])+' ['+qy['unit']+'])'
+                        statunit= '['+qx['unit']+']'
+                    else:
+                        refstr = ('(@ %'+form) % (qx['value'])+' ['+qx['unit']+'])'
+
+            if self._is_all_chan_flagged(i):
+                if not rtnabc:
+                    statvals[i] = None
+                if skip_flaggedrow:
+                    continue
 
             tm = self._gettime(i)
@@ -659 +987 @@
             if len(rows) > 1:
                 # out += ('= %'+form) % (outvec[i]) +'   '+refstr+'\n'
-                out += ('= %'+form) % (statvals[i]) +'   '+refstr+'\n'
+                if statvals[i] is None:
+                    out += ('= None(flagged)') + '   '+refstr+'\n'
+                else:
+                    out += ('= %'+form) % (statvals[i]) +'   '+refstr+'\n'
             else:
                 # out += ('= %'+form) % (outvec[0]) +'   '+refstr+'\n'
-                out += ('= %'+form) % (statvals[0]) +'   '+refstr+'\n'
+                if statvals[0] is None:
+                    out += ('= None(flagged)') + '   '+refstr+'\n'
+                else:
+                    out += ('= %'+form) % (statvals[0]) +'   '+refstr+'\n'
             out +=  sep+"\n"
 
@@ -683 +1017 @@
         asaplog.push(''.join(x), False)
 
+        if skip_flaggedrow:
+            nstatvals = len(statvals)
+            for i in reversed(xrange(nstatvals)):
+                if statvals[i] is None:
+                    del statvals[i]
         return statvals
 
@@ -765 +1104 @@
         return self._get_column( self._gettsysspectrum, row )
 
+    def set_tsys(self, values, row=-1):
+        """\
+        Set the Tsys value(s) of the given 'row' or the whole scantable 
+        (selection).
+        
+        Parameters:
+
+            values:    a scalar or list (if Tsys is a vector) of Tsys value(s)
+            row:       the row number to apply Tsys values to. 
+                       (default all rows)
+                       
+        """
+            
+        if not hasattr(values, "__len__"):
+            values = [values]
+        self._settsys(values, row)
+
     def get_weather(self, row=-1):
         """\
-        Return the weather informations.
+        Return the weather information.
 
         Parameters:
@@ -778 +1134 @@
 
         """
-
+        if row >= len(self):
+            raise IndexError("row out of range")
         values = self._get_column(self._get_weather, row)
         if row > -1:
@@ -788 +1145 @@
             out = []
             for r in values:
-
                 out.append({'temperature': r[0],
                             'pressure': r[1], 'humidity' : r[2],
@@ -1005 +1361 @@
         self._add_history("set_feedtype", vars())
 
+    @asaplog_post_dec
+    def get_doppler(self):
+        """\
+        Get the doppler.
+        """
+        return self._getcoordinfo()[2]
+    
     @asaplog_post_dec
     def set_doppler(self, doppler='RADIO'):
@@ -1471 +1834 @@
 
     @asaplog_post_dec
+    def parse_spw_selection(self, selectstring, restfreq=None, frame=None, doppler=None):
+        """
+        Parse MS type spw/channel selection syntax. 
+
+        Parameters:
+            selectstring : A string expression of spw and channel selection.
+                           Comma-separated expressions mean different spw -
+                           channel combinations. Spws and channel selections
+                           are partitioned by a colon ':'. In a single
+                           selection expression, you can put multiple values
+                           separated by semicolons ';'. Both for spw and
+                           channel selection, allowed cases include single
+                           value, blank('') or asterisk('*') to specify all
+                           available values, two values connected with a
+                           tilde ('~') to specify an inclusive range. Unit
+                           strings for frequency or velocity can be added to
+                           the tilde-connected values. For channel selection
+                           expression, placing a '<' or a '>' is possible to
+                           specify a semi-infinite interval as well.
+
+                     examples:
+                           '' or '*'   = all spws (all channels)
+                           '<2,4~6,9'  = Spws 0,1,4,5,6,9 (all channels)
+                           '3:3~45;60' = channels 3 to 45 and 60 in spw 3
+                           '0~1:2~6,8' = channels 2 to 6 in spws 0,1, and
+                                         all channels in spw8
+                           '1.3~1.5GHz' = all spws whose central frequency
+                                          falls in frequency range between
+                                          1.3GHz and 1.5GHz.
+                           '1.3~1.5GHz:1.3~1.5GHz' = channels which fall
+                                                     between the specified
+                                                     frequency range in spws
+                                                     whose central frequency
+                                                     falls in the specified
+                                                     frequency range.
+                           '1:-200~250km/s' = channels that fall between the
+                                              specified velocity range in
+                                              spw 1.
+            restfreq: the rest frequency.
+                     examples: '115.2712GHz', 115271201800.0
+            frame:   an optional frame type, default 'LSRK'. Valid frames are:
+                     'TOPO', 'LSRD', 'LSRK', 'BARY',
+                     'GEO', 'GALACTO', 'LGROUP', 'CMB'
+            doppler: one of 'RADIO', 'OPTICAL', 'Z', 'BETA', 'GAMMA'
+        Returns:
+        A dictionary of selected (valid) spw and masklist pairs,
+        e.g. {'0': [[50,250],[350,462]], '2': [[100,400],[550,974]]}
+        """
+        if not isinstance(selectstring, str):
+            asaplog.post()
+            asaplog.push("Expression of spw/channel selection must be a string.")
+            asaplog.post("ERROR")
+
+        orig_unit = self.get_unit()
+        self.set_unit('channel')
+        
+        if restfreq is not None:
+            orig_molids = self._getmolidcol_list()
+            set_restfreq(self, restfreq)
+
+        orig_coord = self._getcoordinfo()
+
+        if frame is not None:
+            orig_frame = orig_coord[1]
+            self.set_freqframe(frame)
+
+        if doppler is not None:
+            orig_doppler = orig_coord[2]
+            self.set_doppler(doppler)
+        
+        valid_ifs = self.getifnos()
+
+        comma_sep = selectstring.split(",")
+        res = {}
+
+        for cms_elem in comma_sep:
+            colon_sep = cms_elem.split(":")
+            
+            if (len(colon_sep) > 2):
+                raise RuntimeError("Invalid selection expression: more than two colons!")
+            
+            # parse spw expression and store result in spw_list.
+            # allowed cases include '', '*', 'a', '<a', '>a', 'a~b',
+            # 'a~b*Hz' (where * can be '', 'k', 'M', 'G' etc.),
+            # 'a~b*m/s' (where * can be '' or 'k') and also
+            # several of the above expressions connected with ';'.
+            
+            spw_list = []
+
+            semicolon_sep = colon_sep[0].split(";")
+            
+            for scs_elem in semicolon_sep:
+                scs_elem = scs_elem.strip()
+                
+                lt_sep = scs_elem.split("<")
+                gt_sep = scs_elem.split(">")
+                ti_sep = scs_elem.split("~")
+                
+                lt_sep_length = len(lt_sep)
+                gt_sep_length = len(gt_sep)
+                ti_sep_length = len(ti_sep)
+                
+                len_product = lt_sep_length * gt_sep_length * ti_sep_length
+
+                if (len_product > 2):
+                    # '<', '>' and '~' must not coexist in a single spw expression
+                    
+                    raise RuntimeError("Invalid spw selection.")
+                
+                elif (len_product == 1):
+                    # '', '*', or single spw number.
+                    
+                    if (scs_elem == "") or (scs_elem == "*"):
+                        spw_list = valid_ifs[:] # deep copy
+                    
+                    else: # single number
+                        expr = int(scs_elem)
+                        spw_list.append(expr)
+                        if expr not in valid_ifs:
+                            asaplog.push("Invalid spw given. Ignored.")
+                    
+                else: # (len_product == 2)
+                    # namely, one of '<', '>' or '~' appears just once.
+                    
+                    if (lt_sep_length == 2): # '<a'
+                        if is_number(lt_sep[1]):
+                            no_valid_spw = True
+                            for i in valid_ifs:
+                                if (i < float(lt_sep[1])):
+                                    spw_list.append(i)
+                                    no_valid_spw = False
+
+                            if no_valid_spw:
+                                raise ValueError("Invalid spw selection ('<" + str(lt_sep[1]) + "').")
+                        
+                        else:
+                            raise RuntimeError("Invalid spw selection.")
+                        
+                    elif (gt_sep_length == 2): # '>a'
+                        if is_number(gt_sep[1]):
+                            no_valid_spw = True
+                            for i in valid_ifs:
+                                if (i > float(gt_sep[1])):
+                                    spw_list.append(i)
+                                    no_valid_spw = False
+                            
+                            if no_valid_spw:
+                                raise ValueError("Invalid spw selection ('>" + str(gt_sep[1]) + "').")
+                        
+                        else:
+                            raise RuntimeError("Invalid spw selection.")
+                        
+                    else: # (ti_sep_length == 2) where both boundaries inclusive
+                        expr0 = ti_sep[0].strip()
+                        expr1 = ti_sep[1].strip()
+
+                        if is_number(expr0) and is_number(expr1):
+                            # 'a~b'
+                            expr_pmin = min(float(expr0), float(expr1))
+                            expr_pmax = max(float(expr0), float(expr1))
+                            has_invalid_spw = False
+                            no_valid_spw = True
+                            
+                            for i in valid_ifs:
+                                if (expr_pmin <= i) and (i <= expr_pmax):
+                                    spw_list.append(i)
+                                    no_valid_spw = False
+                                else:
+                                    has_invalid_spw = True
+
+                            if has_invalid_spw:
+                                msg = "Invalid spw is given. Ignored."
+                                asaplog.push(msg)
+                                asaplog.post()
+
+                            if no_valid_spw:
+                                raise ValueError("No valid spw in range ('" + str(expr_pmin) + "~" + str(expr_pmax) + "').")
+                        
+                        elif is_number(expr0) and is_frequency(expr1):
+                            # 'a~b*Hz'
+                            (expr_f0, expr_f1) = get_freq_by_string(expr0, expr1)
+                            expr_fmin = min(expr_f0, expr_f1)
+                            expr_fmax = max(expr_f0, expr_f1)
+                            no_valid_spw = True
+                            
+                            for coord in self._get_coordinate_list():
+                                spw = coord['if']
+                                
+                                """
+                                expr_p0 = coord['coord'].to_pixel(expr_f0)
+                                expr_p1 = coord['coord'].to_pixel(expr_f1)
+                                expr_pmin = min(expr_p0, expr_p1)
+                                expr_pmax = max(expr_p0, expr_p1)
+                                
+                                pmin = 0.0
+                                pmax = float(self.nchan(spw) - 1)
+
+                                if ((expr_pmax - pmin)*(expr_pmin - pmax) <= 0.0):
+                                    spw_list.append(spw)
+                                    no_valid_spw = False
+                                """
+
+                                crd = coord['coord']
+                                fhead = crd.to_frequency(0)
+                                ftail = crd.to_frequency(self.nchan(spw) - 1)
+                                fcen  = (fhead + ftail) / 2.0
+
+                                if ((expr_fmin <= fcen) and (fcen <= expr_fmax)):
+                                    spw_list.append(spw)
+                                    no_valid_spw = False
+                                
+                            if no_valid_spw:
+                                raise ValueError("No valid spw in range ('" + str(expr0) + "~" + str(expr1) + "').")
+                                
+                        elif is_number(expr0) and is_velocity(expr1):
+                            # 'a~b*m/s'
+                            (expr_v0, expr_v1) = get_velocity_by_string(expr0, expr1)
+                            expr_vmin = min(expr_v0, expr_v1)
+                            expr_vmax = max(expr_v0, expr_v1)
+                            no_valid_spw = True
+                            
+                            for coord in self._get_coordinate_list():
+                                spw = coord['if']
+
+                                """
+                                pmin = 0.0
+                                pmax = float(self.nchan(spw) - 1)
+                                
+                                vel0 = coord['coord'].to_velocity(pmin)
+                                vel1 = coord['coord'].to_velocity(pmax)
+                                
+                                vmin = min(vel0, vel1)
+                                vmax = max(vel0, vel1)
+
+                                if ((expr_vmax - vmin)*(expr_vmin - vmax) <= 0.0):
+                                    spw_list.append(spw)
+                                    no_valid_spw = False
+                                """
+
+                                crd = coord['coord']
+                                fhead = crd.to_frequency(0)
+                                ftail = crd.to_frequency(self.nchan(spw) - 1)
+                                fcen  = (fhead + ftail) / 2.0
+                                vcen  = crd.to_velocity(crd.to_pixel(fcen))
+
+                                if ((expr_vmin <= vcen) and (vcen <= expr_vmax)):
+                                    spw_list.append(spw)
+                                    no_valid_spw = False
+
+                            if no_valid_spw:
+                                raise ValueError("No valid spw in range ('" + str(expr0) + "~" + str(expr1) + "').")
+                            
+                        else:
+                            # cases such as 'aGHz~bkm/s' are not allowed now
+                            raise RuntimeError("Invalid spw selection.")
+
+            # check spw list and remove invalid ones.
+            # if no valid spw left, emit ValueError.
+            if len(spw_list) == 0:
+                raise ValueError("No valid spw in given range.")
+            
+            # parse channel expression and store the result in crange_list.
+            # allowed cases include '', 'a~b', 'a*Hz~b*Hz' (where * can be
+            # '', 'k', 'M', 'G' etc.), 'a*m/s~b*m/s' (where * can be '' or 'k')
+            # and also several of the above expressions connected with ';'.
+            
+            for spw in spw_list:
+                pmin = 0.0
+                pmax = float(self.nchan(spw) - 1)
+
+                molid = self._getmolidcol_list()[self.get_first_rowno_by_if(spw)]
+                
+                if (len(colon_sep) == 1):
+                    # no expression for channel selection, 
+                    # which means all channels are to be selected.
+                    crange_list = [[pmin, pmax]]
+                
+                else: # (len(colon_sep) == 2)
+                    crange_list = []
+                    
+                    found = False
+                    for i in self._get_coordinate_list():
+                        if (i['if'] == spw):
+                            coord = i['coord']
+                            found = True
+                            break
+
+                    if found:
+                        semicolon_sep = colon_sep[1].split(";")
+                        for scs_elem in semicolon_sep:
+                            scs_elem = scs_elem.strip()
+
+                            ti_sep = scs_elem.split("~")
+                            ti_sep_length = len(ti_sep)
+
+                            if (ti_sep_length > 2):
+                                raise RuntimeError("Invalid channel selection.")
+                        
+                            elif (ti_sep_length == 1):
+                                if (scs_elem == "") or (scs_elem == "*"):
+                                    # '' and '*' for all channels
+                                    crange_list = [[pmin, pmax]]
+                                    break
+                                elif (is_number(scs_elem)):
+                                    # single channel given
+                                    crange_list.append([float(scs_elem), float(scs_elem)])
+                                else:
+                                    raise RuntimeError("Invalid channel selection.")
+
+                            else: #(ti_sep_length == 2)
+                                expr0 = ti_sep[0].strip()
+                                expr1 = ti_sep[1].strip()
+
+                                if is_number(expr0) and is_number(expr1):
+                                    # 'a~b'
+                                    expr_pmin = min(float(expr0), float(expr1))
+                                    expr_pmax = max(float(expr0), float(expr1))
+
+                                elif is_number(expr0) and is_frequency(expr1):
+                                    # 'a~b*Hz'
+                                    (expr_f0, expr_f1) = get_freq_by_string(expr0, expr1)
+                                    expr_p0 = coord.to_pixel(expr_f0)
+                                    expr_p1 = coord.to_pixel(expr_f1)
+                                    expr_pmin = min(expr_p0, expr_p1)
+                                    expr_pmax = max(expr_p0, expr_p1)
+
+                                elif is_number(expr0) and is_velocity(expr1):
+                                    # 'a~b*m/s'
+                                    restf = self.get_restfreqs()[molid][0]
+                                    (expr_v0, expr_v1) = get_velocity_by_string(expr0, expr1)
+                                    dppl = self.get_doppler()
+                                    expr_f0 = get_frequency_by_velocity(restf, expr_v0, dppl)
+                                    expr_f1 = get_frequency_by_velocity(restf, expr_v1, dppl)
+                                    expr_p0 = coord.to_pixel(expr_f0)
+                                    expr_p1 = coord.to_pixel(expr_f1)
+                                    expr_pmin = min(expr_p0, expr_p1)
+                                    expr_pmax = max(expr_p0, expr_p1)
+                            
+                                else:
+                                    # cases such as 'aGHz~bkm/s' are not allowed now
+                                    raise RuntimeError("Invalid channel selection.")
+
+                                cmin = max(pmin, expr_pmin)
+                                cmax = min(pmax, expr_pmax)
+                                # if the given range of channel selection has overwrap with
+                                # that of current spw, output the overwrap area.
+                                if (cmin <= cmax):
+                                    cmin = float(int(cmin + 0.5))
+                                    cmax = float(int(cmax + 0.5))
+                                    crange_list.append([cmin, cmax])
+
+                    if (len(crange_list) == 0):
+                        crange_list.append([])
+
+                if (len(crange_list[0]) > 0):
+                    if res.has_key(spw):
+                        res[spw].extend(crange_list)
+                    else:
+                        res[spw] = crange_list
+
+        for spw in res.keys():
+            if spw in valid_ifs:
+                # remove duplicated channel ranges
+                for i in reversed(xrange(len(res[spw]))):
+                    for j in xrange(i):
+                        if ((res[spw][i][0]-res[spw][j][1])*(res[spw][i][1]-res[spw][j][0]) <= 0) or \
+                            (min(abs(res[spw][i][0]-res[spw][j][1]),abs(res[spw][j][0]-res[spw][i][1])) == 1):
+                            asaplog.post()
+                            merge_warn_mesg = "Spw " + str(spw) + ": overwrapping channel ranges are merged."
+                            asaplog.push(merge_warn_mesg)
+                            asaplog.post('WARN')
+                            res[spw][j][0] = min(res[spw][i][0], res[spw][j][0])
+                            res[spw][j][1] = max(res[spw][i][1], res[spw][j][1])
+                            res[spw].pop(i)
+                            break
+            else:
+                del res[spw]
+
+        if len(res) == 0:
+            raise RuntimeError("No valid spw.")
+        
+        # restore original values
+        self.set_unit(orig_unit)
+        if restfreq is not None:
+            self._setmolidcol_list(orig_molids)
+        if frame is not None:
+            self.set_freqframe(orig_frame)
+        if doppler is not None:
+            self.set_doppler(orig_doppler)
+        
+        return res
+    
+    @asaplog_post_dec
+    def get_first_rowno_by_if(self, ifno):
+        found = False
+        for irow in xrange(self.nrow()):
+            if (self.getif(irow) == ifno):
+                res = irow
+                found = True
+                break
+
+        if not found: raise RuntimeError("No valid spw.")
+        
+        return res
+
+    @asaplog_post_dec
+    def _get_coordinate_list(self):
+        res = []
+        spws = self.getifnos()
+        for spw in spws:
+            elem = {}
+            elem['if']    = spw
+            elem['coord'] = self.get_coordinate(self.get_first_rowno_by_if(spw))
+            res.append(elem)
+
+        return res
+    
+    @asaplog_post_dec
     def parse_maskexpr(self, maskstring):
         """
@@ -1501 +2282 @@
             maskstring = str(valid_ifs)[1:-1]
         ## split each selection "IF range[:CHAN range]"
-        sellist = maskstring.split(',')
+        # split maskstring by "<spaces>,<spaces>"
+        comma_sep = re.compile('\s*,\s*')
+        sellist = comma_sep.split(maskstring)
+        # separator by "<spaces>:<spaces>"
+        collon_sep = re.compile('\s*:\s*')
         for currselstr in sellist:
-            selset = currselstr.split(':')
+            selset = collon_sep.split(currselstr)
             # spw and mask string (may include ~, < or >)
             spwmasklist = self._parse_selection(selset[0], typestr='integer',
@@ -1637 +2422 @@
         if smode == 'r':
             minidx = 0
-            maxidx = self.nrow()
+            maxidx = self.nrow()-1
         else:
             idx = getattr(self,"get"+mode+"nos")()
@@ -1643 +2428 @@
             maxidx = max(idx)
             del idx
-        sellist = selexpr.split(',')
+        # split selexpr by "<spaces>,<spaces>"
+        comma_sep = re.compile('\s*,\s*')
+        sellist = comma_sep.split(selexpr)
         idxlist = []
         for currselstr in sellist:
@@ -1651 +2438 @@
             for thelist in currlist:
                 idxlist += range(thelist[0],thelist[1]+1)
+        # remove duplicated elements after first ones
+        for i in reversed(xrange(len(idxlist))):
+            if idxlist.index(idxlist[i]) < i:
+                idxlist.pop(i)
+
+        # remove elements outside range [minidx, maxidx] for smode='r'
+        if smode == 'r':
+            for i in reversed(xrange(len(idxlist))):
+                if (idxlist[i] < minidx) or (idxlist[i] > maxidx):
+                    idxlist.pop(i)
+        
         msg = "Selected %s: %s" % (mode.upper()+"NO", str(idxlist))
         asaplog.push(msg)
@@ -1676 +2474 @@
             --> returns [[0.,2.5],[5.0,7.0],[9.,9.]]
         """
-        selgroups = selstr.split(';')
+        # split selstr by '<spaces>;<spaces>'
+        semi_sep = re.compile('\s*;\s*')
+        selgroups = semi_sep.split(selstr)
         sellists = []
         if typestr.lower().startswith('int'):
@@ -1685 +2485 @@
         
         for currsel in  selgroups:
+            if currsel.strip() == '*' or len(currsel.strip()) == 0:
+                minsel = minval
+                maxsel = maxval
             if currsel.find('~') > 0:
                 # val0 <= x <= val1
                 minsel = formatfunc(currsel.split('~')[0].strip())
-                maxsel = formatfunc(currsel.split('~')[1].strip()) 
+                maxsel = formatfunc(currsel.split('~')[1].strip())
             elif currsel.strip().find('<=') > -1:
                 bound = currsel.split('<=')
@@ -1907 +2710 @@
 
     @asaplog_post_dec
-    def history(self, filename=None):
+    def history(self, filename=None, nrows=-1, start=0):
         """\
         Print the history. Optionally to a file.
@@ -1916 +2719 @@
 
         """
-        hist = list(self._gethistory())
+        n = self._historylength()
+        if nrows == -1:
+            nrows = n
+        if start+nrows > n:
+            nrows = nrows-start
+        if n > 1000 and nrows == n:
+            nrows = 1000
+            start = n-1000
+            asaplog.push("Warning: History has {0} entries. Displaying last "
+                         "1000".format(n))
+        hist = list(self._gethistory(nrows, start))
         out = "-"*80
         for h in hist:
-            if h.startswith("---"):
-                out = "\n".join([out, h])
+            if not h.strip():
+                continue
+            if h.find("---") >-1:
+                continue
             else:
                 items = h.split("##")
@@ -1933 +2748 @@
                     s = i.split("=")
                     out += "\n   %s = %s" % (s[0], s[1])
-                out = "\n".join([out, "-"*80])
+                out = "\n".join([out, "*"*80])
         if filename is not None:
             if filename is "":
                 filename = 'scantable_history.txt'
-            import os
             filename = os.path.expandvars(os.path.expanduser(filename))
             if not os.path.isdir(filename):
@@ -1952 +2766 @@
     #
     @asaplog_post_dec
-    def average_time(self, mask=None, scanav=False, weight='tint', align=False):
+    def average_time(self, mask=None, scanav=False, weight='tint', align=False,
+                     avmode="NONE"):
         """\
         Return the (time) weighted average of a scan. Scans will be averaged
@@ -1980 +2795 @@
             align:    align the spectra in velocity before averaging. It takes
                       the time of the first spectrum as reference time.
+            avmode:   'SOURCE' - also select by source name -  or 
+                      'NONE' (default). Not applicable for scanav=True or
+                      weight=median
 
         Example::
@@ -1990 +2808 @@
         weight = weight or 'TINT'
         mask = mask or ()
-        scanav = (scanav and 'SCAN') or 'NONE'
+        scanav = (scanav and 'SCAN') or avmode.upper()
         scan = (self, )
 
         if align:
             scan = (self.freq_align(insitu=False), )
+            asaplog.push("Note: Alignment is don on a source-by-source basis")
+            asaplog.push("Note: Averaging (by default) is not")
+            # we need to set it to SOURCE averaging here            
         s = None
         if weight.upper() == 'MEDIAN':
@@ -2539 +3360 @@
             raise RuntimeError(msg)
 
-
-    @asaplog_post_dec
-    def sinusoid_baseline(self, insitu=None, mask=None, applyfft=None, 
+    @asaplog_post_dec
+    def apply_bltable(self, insitu=None, retfitres=None, inbltable=None, outbltable=None, overwrite=None):
+        """\
+        Subtract baseline based on parameters written in Baseline Table.
+
+        Parameters:
+            insitu:        if True, baseline fitting/subtraction is done
+                           in-situ. If False, a new scantable with
+                           baseline subtracted is returned. Actually,
+                           format of the returned value depends on both 
+                           insitu and retfitres (see below). 
+                           The default is taken from .asaprc (False)
+            retfitres:     if True, the results of baseline fitting (i.e.,
+                           coefficients and rms) are returned.
+                           default is False.
+                           The format of the returned value of this
+                           function varies as follows:
+                           (1) in case insitu=True and retfitres=True:
+                               fitting result.
+                           (2) in case insitu=True and retfitres=False:
+                               None.
+                           (3) in case insitu=False and retfitres=True:
+                               a dictionary containing a new scantable
+                               (with baseline subtracted) and the fitting
+                               results.
+                           (4) in case insitu=False and retfitres=False:
+                               a new scantable (with baseline subtracted).
+            inbltable:     name of input baseline table. The row number of
+                           scantable and that of inbltable must be
+                           identical.
+            outbltable:    name of output baseline table where baseline 
+                           parameters and fitting results recorded.
+                           default is ''(no output).
+            overwrite:     if True when an existing baseline table is
+                           specified for outbltable, overwrites it. 
+                           Otherwise there is no harm.
+                           default is False.
+        """
+
+        try:
+            varlist = vars()
+            if retfitres      is None: retfitres      = False
+            if inbltable      is None: raise ValueError("bltable missing.")
+            if outbltable     is None: outbltable     = ''
+            if overwrite      is None: overwrite      = False
+
+            if insitu is None: insitu = rcParams['insitu']
+            if insitu:
+                workscan = self
+            else:
+                workscan = self.copy()
+
+            sres = workscan._apply_bltable(inbltable,
+                                           retfitres,
+                                           outbltable,
+                                           os.path.exists(outbltable),
+                                           overwrite)
+            if retfitres: res = parse_fitresult(sres)
+
+            workscan._add_history('apply_bltable', varlist)
+
+            if insitu:
+                self._assign(workscan)
+                if retfitres:
+                    return res
+                else:
+                    return None
+            else:
+                if retfitres:
+                    return {'scantable': workscan, 'fitresults': res}
+                else:
+                    return workscan
+        
+        except RuntimeError, e:
+            raise_fitting_failure_exception(e)
+
+    @asaplog_post_dec
+    def sub_baseline(self, insitu=None, retfitres=None, blinfo=None, bltable=None, overwrite=None):
+        """\
+        Subtract baseline based on parameters written in the input list. 
+
+        Parameters:
+            insitu:        if True, baseline fitting/subtraction is done
+                           in-situ. If False, a new scantable with
+                           baseline subtracted is returned. Actually,
+                           format of the returned value depends on both 
+                           insitu and retfitres (see below). 
+                           The default is taken from .asaprc (False)
+            retfitres:     if True, the results of baseline fitting (i.e.,
+                           coefficients and rms) are returned.
+                           default is False.
+                           The format of the returned value of this
+                           function varies as follows:
+                           (1) in case insitu=True and retfitres=True:
+                               fitting result.
+                           (2) in case insitu=True and retfitres=False:
+                               None.
+                           (3) in case insitu=False and retfitres=True:
+                               a dictionary containing a new scantable
+                               (with baseline subtracted) and the fitting
+                               results.
+                           (4) in case insitu=False and retfitres=False:
+                               a new scantable (with baseline subtracted).
+            blinfo:        baseline parameter set stored in a dictionary
+                           or a list of dictionary. Each dictionary 
+                           corresponds to each spectrum and must contain
+                           the following keys and values:
+                             'row': row number, 
+                             'blfunc': function name. available ones include
+                                       'poly', 'chebyshev', 'cspline' and
+                                       'sinusoid', 
+                             'order': maximum order of polynomial. needed
+                                      if blfunc='poly' or 'chebyshev', 
+                             'npiece': number or piecewise polynomial.
+                                       needed if blfunc='cspline', 
+                             'nwave': a list of sinusoidal wave numbers.
+                                      needed if blfunc='sinusoid', and 
+                             'masklist': min-max windows for channel mask.
+                                         the specified ranges will be used
+                                         for fitting. 
+            bltable:       name of output baseline table where baseline
+                           parameters and fitting results recorded.
+                           default is ''(no output).
+            overwrite:     if True when an existing baseline table is
+                           specified for bltable, overwrites it. 
+                           Otherwise there is no harm.
+                           default is False.
+                           
+        Example:
+            sub_baseline(blinfo=[{'row':0, 'blfunc':'poly', 'order':5,
+                                  'masklist':[[10,350],[352,510]]},
+                                 {'row':1, 'blfunc':'cspline', 'npiece':3,
+                                  'masklist':[[3,16],[19,404],[407,511]]}
+                                  ])
+
+                the first spectrum (row=0) will be fitted with polynomial
+                of order=5 and the next one (row=1) will be fitted with cubic
+                spline consisting of 3 pieces. 
+        """
+
+        try:
+            varlist = vars()
+            if retfitres      is None: retfitres      = False
+            if blinfo         is None: blinfo         = []
+            if bltable        is None: bltable        = ''
+            if overwrite      is None: overwrite      = False
+
+            if insitu is None: insitu = rcParams['insitu']
+            if insitu:
+                workscan = self
+            else:
+                workscan = self.copy()
+
+            nrow = workscan.nrow()
+
+            in_blinfo = pack_blinfo(blinfo=blinfo, maxirow=nrow)
+
+            sres = workscan._sub_baseline(in_blinfo,
+                                          retfitres,
+                                          bltable,
+                                          os.path.exists(bltable),
+                                          overwrite)
+            if retfitres: res = parse_fitresult(sres)
+            
+            workscan._add_history('sub_baseline', varlist)
+
+            if insitu:
+                self._assign(workscan)
+                if retfitres:
+                    return res
+                else:
+                    return None
+            else:
+                if retfitres:
+                    return {'scantable': workscan, 'fitresults': res}
+                else:
+                    return workscan
+
+        except RuntimeError, e:
+            raise_fitting_failure_exception(e)
+
+    @asaplog_post_dec
+    def calc_aic(self, value=None, blfunc=None, order=None, mask=None,
+                 whichrow=None, uselinefinder=None, edge=None,
+                 threshold=None, chan_avg_limit=None):
+        """\
+        Calculates and returns model selection criteria for a specified
+        baseline model and a given spectrum data. 
+        Available values include Akaike Information Criterion (AIC), the
+        corrected Akaike Information Criterion (AICc) by Sugiura(1978),
+        Bayesian Information Criterion (BIC) and the Generalised Cross
+        Validation (GCV). 
+
+        Parameters:
+            value:         name of model selection criteria to calculate.
+                           available ones include 'aic', 'aicc', 'bic' and
+                           'gcv'. default is 'aicc'.
+            blfunc:        baseline function name. available ones include 
+                           'chebyshev', 'cspline' and 'sinusoid'.
+                           default is 'chebyshev'.
+            order:         parameter for basline function. actually stands for
+                           order of polynomial (order) for 'chebyshev',
+                           number of spline pieces (npiece) for 'cspline' and
+                           maximum wave number for 'sinusoid', respectively.
+                           default is 5 (which is also the default order value
+                           for [auto_]chebyshev_baseline()).
+            mask:          an optional mask. default is [].
+            whichrow:      row number. default is 0 (the first row)
+            uselinefinder: use sd.linefinder() to flag out line regions
+                           default is True.
+            edge:           an optional number of channel to drop at
+                            the edge of spectrum. If only one value is
+                            specified, the same number will be dropped
+                            from both sides of the spectrum. Default
+                            is to keep all channels. Nested tuples
+                            represent individual edge selection for
+                            different IFs (a number of spectral channels
+                            can be different)
+                            default is (0, 0).
+            threshold:      the threshold used by line finder. It is
+                            better to keep it large as only strong lines
+                            affect the baseline solution.
+                            default is 3.
+            chan_avg_limit: a maximum number of consequtive spectral
+                            channels to average during the search of
+                            weak and broad lines. The default is no
+                            averaging (and no search for weak lines).
+                            If such lines can affect the fitted baseline
+                            (e.g. a high order polynomial is fitted),
+                            increase this parameter (usually values up
+                            to 8 are reasonable). Most users of this
+                            method should find the default value sufficient.
+                            default is 1.
+
+        Example:
+            aic = scan.calc_aic(blfunc='chebyshev', order=5, whichrow=0)
+        """
+
+        try:
+            varlist = vars()
+
+            if value          is None: value          = 'aicc'
+            if blfunc         is None: blfunc         = 'chebyshev'
+            if order          is None: order          = 5
+            if mask           is None: mask           = []
+            if whichrow       is None: whichrow       = 0
+            if uselinefinder  is None: uselinefinder  = True
+            if edge           is None: edge           = (0, 0)
+            if threshold      is None: threshold      = 3
+            if chan_avg_limit is None: chan_avg_limit = 1
+
+            return self._calc_aic(value, blfunc, order, mask,
+                                  whichrow, uselinefinder, edge,
+                                  threshold, chan_avg_limit)
+            
+        except RuntimeError, e:
+            raise_fitting_failure_exception(e)
+
+    @asaplog_post_dec
+    def sinusoid_baseline(self, mask=None, applyfft=None, 
                           fftmethod=None, fftthresh=None,
-                          addwn=None, rejwn=None, clipthresh=None,
-                          clipniter=None, plot=None,
-                          getresidual=None, showprogress=None, 
-                          minnrow=None, outlog=None, blfile=None, csvformat=None):
+                          addwn=None, rejwn=None,
+                          insitu=None,
+                          clipthresh=None, clipniter=None,
+                          plot=None,
+                          getresidual=None,
+                          showprogress=None, minnrow=None,
+                          outlog=None,
+                          blfile=None, csvformat=None,
+                          bltable=None):
         """\
         Return a scan which has been baselined (all rows) with sinusoidal 
@@ -2552 +3635 @@
 
         Parameters:
-            insitu:        if False a new scantable is returned.
-                           Otherwise, the scaling is done in-situ
-                           The default is taken from .asaprc (False)
             mask:          an optional mask
             applyfft:      if True use some method, such as FFT, to find
@@ -2586 +3666 @@
                            wave numbers which are specified both in addwn
                            and rejwn will NOT be used. default is []. 
+            insitu:        if False a new scantable is returned.
+                           Otherwise, the scaling is done in-situ
+                           The default is taken from .asaprc (False)
             clipthresh:    Clipping threshold. (default is 3.0, unit: sigma)
             clipniter:     maximum number of iteration of 'clipthresh'-sigma 
@@ -2605 +3688 @@
                            (default is '': no file/logger output)
             csvformat:     if True blfile is csv-formatted, default is False.
+            bltable:       name of a baseline table where fitting results
+                           (coefficients, rms, etc.) are to be written.
+                           if given, fitting results will NOT be output to
+                           scantable (insitu=True) or None will be 
+                           returned (insitu=False). 
+                           (default is "": no table output)
 
         Example:
@@ -2626 +3715 @@
                 workscan = self.copy()
             
-            #if mask          is None: mask          = [True for i in xrange(workscan.nchan())]
             if mask          is None: mask          = []
             if applyfft      is None: applyfft      = True
@@ -2642 +3730 @@
             if blfile        is None: blfile        = ''
             if csvformat     is None: csvformat     = False
-
-            if csvformat:
-                scsvformat = "T"
-            else:
-                scsvformat = "F"
+            if bltable       is None: bltable       = ''
+
+            sapplyfft = 'true' if applyfft else 'false'
+            fftinfo = ','.join([sapplyfft, fftmethod.lower(), str(fftthresh).lower()])
+
+            scsvformat = 'T' if csvformat else 'F'
 
             #CURRENTLY, PLOT=true is UNAVAILABLE UNTIL sinusoidal fitting is implemented as a fitter method. 
-            workscan._sinusoid_baseline(mask, applyfft, fftmethod.lower(), 
-                                        str(fftthresh).lower(), 
+            workscan._sinusoid_baseline(mask,
+                                        fftinfo, 
+                                        #applyfft, fftmethod.lower(), 
+                                        #str(fftthresh).lower(), 
                                         workscan._parse_wn(addwn), 
                                         workscan._parse_wn(rejwn),
@@ -2657 +3748 @@
                                         pack_progress_params(showprogress, 
                                                              minnrow),
-                                        outlog, scsvformat+blfile)
+                                        outlog, scsvformat+blfile,
+                                        bltable)
             workscan._add_history('sinusoid_baseline', varlist)
-            
-            if insitu:
-                self._assign(workscan)
+
+            if bltable == '':
+                if insitu:
+                    self._assign(workscan)
+                else:
+                    return workscan
             else:
-                return workscan
+                if not insitu:
+                    return None
             
         except RuntimeError, e:
@@ -2670 +3766 @@
 
     @asaplog_post_dec
-    def auto_sinusoid_baseline(self, insitu=None, mask=None, applyfft=None, 
+    def auto_sinusoid_baseline(self, mask=None, applyfft=None, 
                                fftmethod=None, fftthresh=None,
-                               addwn=None, rejwn=None, clipthresh=None,
-                               clipniter=None, edge=None, threshold=None,
-                               chan_avg_limit=None, plot=None,
-                               getresidual=None, showprogress=None, 
-                               minnrow=None, outlog=None, blfile=None, csvformat=None):
+                               addwn=None, rejwn=None,
+                               insitu=None,
+                               clipthresh=None, clipniter=None,
+                               edge=None, threshold=None, chan_avg_limit=None,
+                               plot=None,
+                               getresidual=None,
+                               showprogress=None, minnrow=None,
+                               outlog=None,
+                               blfile=None, csvformat=None,
+                               bltable=None):
         """\
         Return a scan which has been baselined (all rows) with sinusoidal 
@@ -2684 +3785 @@
 
         Parameters:
-            insitu:         if False a new scantable is returned.
-                            Otherwise, the scaling is done in-situ
-                            The default is taken from .asaprc (False)
             mask:           an optional mask retreived from scantable
             applyfft:       if True use some method, such as FFT, to find
@@ -2718 +3816 @@
                             wave numbers which are specified both in addwn
                             and rejwn will NOT be used. default is []. 
+            insitu:         if False a new scantable is returned.
+                            Otherwise, the scaling is done in-situ
+                            The default is taken from .asaprc (False)
             clipthresh:     Clipping threshold. (default is 3.0, unit: sigma)
             clipniter:      maximum number of iteration of 'clipthresh'-sigma 
@@ -2757 +3858 @@
                             (default is "": no file/logger output)
             csvformat:      if True blfile is csv-formatted, default is False.
+            bltable:        name of a baseline table where fitting results
+                            (coefficients, rms, etc.) are to be written.
+                            if given, fitting results will NOT be output to
+                            scantable (insitu=True) or None will be 
+                            returned (insitu=False). 
+                            (default is "": no table output)
 
         Example:
@@ -2775 +3882 @@
                 workscan = self.copy()
             
-            #if mask           is None: mask           = [True for i in xrange(workscan.nchan())]
             if mask           is None: mask           = []
             if applyfft       is None: applyfft       = True
@@ -2794 +3900 @@
             if blfile         is None: blfile         = ''
             if csvformat      is None: csvformat      = False
-
-            if csvformat:
-                scsvformat = "T"
-            else:
-                scsvformat = "F"
+            if bltable        is None: bltable        = ''
+
+            sapplyfft = 'true' if applyfft else 'false'
+            fftinfo = ','.join([sapplyfft, fftmethod.lower(), str(fftthresh).lower()])
+
+            scsvformat = 'T' if csvformat else 'F'
 
             #CURRENTLY, PLOT=true is UNAVAILABLE UNTIL sinusoidal fitting is implemented as a fitter method. 
-            workscan._auto_sinusoid_baseline(mask, applyfft, 
-                                             fftmethod.lower(), 
-                                             str(fftthresh).lower(), 
+            workscan._auto_sinusoid_baseline(mask,
+                                             fftinfo, 
                                              workscan._parse_wn(addwn), 
                                              workscan._parse_wn(rejwn), 
@@ -2812 +3918 @@
                                              pack_progress_params(showprogress,
                                                                   minnrow), 
-                                             outlog, scsvformat+blfile)
+                                             outlog, scsvformat+blfile, bltable)
             workscan._add_history("auto_sinusoid_baseline", varlist)
-            
-            if insitu:
-                self._assign(workscan)
+
+            if bltable == '':
+                if insitu:
+                    self._assign(workscan)
+                else:
+                    return workscan
             else:
-                return workscan
+                if not insitu:
+                    return None
             
         except RuntimeError, e:
@@ -2824 +3934 @@
 
     @asaplog_post_dec
-    def cspline_baseline(self, insitu=None, mask=None, npiece=None, 
+    def cspline_baseline(self, mask=None, npiece=None, insitu=None, 
                          clipthresh=None, clipniter=None, plot=None, 
                          getresidual=None, showprogress=None, minnrow=None, 
-                         outlog=None, blfile=None, csvformat=None):
+                         outlog=None, blfile=None, csvformat=None,
+                         bltable=None):
         """\
         Return a scan which has been baselined (all rows) by cubic spline 
@@ -2833 +3944 @@
 
         Parameters:
+            mask:         An optional mask
+            npiece:       Number of pieces. (default is 2)
             insitu:       If False a new scantable is returned.
                           Otherwise, the scaling is done in-situ
                           The default is taken from .asaprc (False)
-            mask:         An optional mask
-            npiece:       Number of pieces. (default is 2)
             clipthresh:   Clipping threshold. (default is 3.0, unit: sigma)
             clipniter:    maximum number of iteration of 'clipthresh'-sigma 
@@ -2857 +3968 @@
                           (default is "": no file/logger output)
             csvformat:    if True blfile is csv-formatted, default is False.
+            bltable:      name of a baseline table where fitting results
+                          (coefficients, rms, etc.) are to be written.
+                          if given, fitting results will NOT be output to
+                          scantable (insitu=True) or None will be 
+                          returned (insitu=False). 
+                          (default is "": no table output)
 
         Example:
@@ -2878 +3995 @@
                 workscan = self.copy()
 
-            #if mask         is None: mask         = [True for i in xrange(workscan.nchan())]
             if mask         is None: mask         = []
             if npiece       is None: npiece       = 2
@@ -2889 +4005 @@
             if outlog       is None: outlog       = False
             if blfile       is None: blfile       = ''
-            if csvformat     is None: csvformat     = False
-
-            if csvformat:
-                scsvformat = "T"
-            else:
-                scsvformat = "F"
+            if csvformat    is None: csvformat    = False
+            if bltable      is None: bltable      = ''
+
+            scsvformat = 'T' if csvformat else 'F'
 
             #CURRENTLY, PLOT=true UNAVAILABLE UNTIL cubic spline fitting is implemented as a fitter method. 
-            workscan._cspline_baseline(mask, npiece, clipthresh, clipniter, 
+            workscan._cspline_baseline(mask, npiece,
+                                       clipthresh, clipniter, 
                                        getresidual, 
                                        pack_progress_params(showprogress, 
                                                             minnrow),
-                                       outlog, scsvformat+blfile)
+                                       outlog, scsvformat+blfile,
+                                       bltable)
             workscan._add_history("cspline_baseline", varlist)
-            
-            if insitu:
-                self._assign(workscan)
+
+            if bltable == '':
+                if insitu:
+                    self._assign(workscan)
+                else:
+                    return workscan
             else:
-                return workscan
+                if not insitu:
+                    return None
             
         except RuntimeError, e:
@@ -2913 +4033 @@
 
     @asaplog_post_dec
-    def auto_cspline_baseline(self, insitu=None, mask=None, npiece=None, 
+    def auto_cspline_baseline(self, mask=None, npiece=None, insitu=None, 
                               clipthresh=None, clipniter=None,
                               edge=None, threshold=None, chan_avg_limit=None, 
                               getresidual=None, plot=None,
                               showprogress=None, minnrow=None, outlog=None,
-                              blfile=None, csvformat=None):
+                              blfile=None, csvformat=None, bltable=None):
         """\
         Return a scan which has been baselined (all rows) by cubic spline
@@ -2926 +4046 @@
 
         Parameters:
+            mask:           an optional mask retreived from scantable
+            npiece:         Number of pieces. (default is 2)
             insitu:         if False a new scantable is returned.
                             Otherwise, the scaling is done in-situ
                             The default is taken from .asaprc (False)
-            mask:           an optional mask retreived from scantable
-            npiece:         Number of pieces. (default is 2)
             clipthresh:     Clipping threshold. (default is 3.0, unit: sigma)
             clipniter:      maximum number of iteration of 'clipthresh'-sigma 
@@ -2970 +4090 @@
                             (default is "": no file/logger output)
             csvformat:      if True blfile is csv-formatted, default is False.
+            bltable:        name of a baseline table where fitting results
+                            (coefficients, rms, etc.) are to be written.
+                            if given, fitting results will NOT be output to
+                            scantable (insitu=True) or None will be 
+                            returned (insitu=False). 
+                            (default is "": no table output)
 
         Example:
@@ -3003 +4129 @@
             if blfile         is None: blfile         = ''
             if csvformat      is None: csvformat      = False
-
-            if csvformat:
-                scsvformat = "T"
-            else:
-                scsvformat = "F"
+            if bltable        is None: bltable        = ''
+
+            scsvformat = 'T' if csvformat else 'F'
 
             #CURRENTLY, PLOT=true UNAVAILABLE UNTIL cubic spline fitting is implemented as a fitter method.
-            workscan._auto_cspline_baseline(mask, npiece, clipthresh, 
-                                            clipniter, 
+            workscan._auto_cspline_baseline(mask, npiece,
+                                            clipthresh, clipniter, 
                                             normalise_edge_param(edge), 
                                             threshold, 
@@ -3017 +4141 @@
                                             pack_progress_params(showprogress, 
                                                                  minnrow), 
-                                            outlog, scsvformat+blfile)
+                                            outlog,
+                                            scsvformat+blfile,
+                                            bltable)
             workscan._add_history("auto_cspline_baseline", varlist)
-            
-            if insitu:
-                self._assign(workscan)
+
+            if bltable == '':
+                if insitu:
+                    self._assign(workscan)
+                else:
+                    return workscan
             else:
-                return workscan
+                if not insitu:
+                    return None
             
         except RuntimeError, e:
@@ -3029 +4159 @@
 
     @asaplog_post_dec
-    def chebyshev_baseline(self, insitu=None, mask=None, order=None, 
+    def chebyshev_baseline(self, mask=None, order=None, insitu=None, 
                            clipthresh=None, clipniter=None, plot=None, 
                            getresidual=None, showprogress=None, minnrow=None, 
-                           outlog=None, blfile=None, csvformat=None):
+                           outlog=None, blfile=None, csvformat=None,
+                           bltable=None):
         """\
         Return a scan which has been baselined (all rows) by Chebyshev polynomials.
 
         Parameters:
+            mask:         An optional mask
+            order:        the maximum order of Chebyshev polynomial (default is 5)
             insitu:       If False a new scantable is returned.
                           Otherwise, the scaling is done in-situ
                           The default is taken from .asaprc (False)
-            mask:         An optional mask
-            order:        the maximum order of Chebyshev polynomial (default is 5)
             clipthresh:   Clipping threshold. (default is 3.0, unit: sigma)
             clipniter:    maximum number of iteration of 'clipthresh'-sigma 
@@ -3061 +4192 @@
                           (default is "": no file/logger output)
             csvformat:    if True blfile is csv-formatted, default is False.
+            bltable:      name of a baseline table where fitting results
+                          (coefficients, rms, etc.) are to be written.
+                          if given, fitting results will NOT be output to
+                          scantable (insitu=True) or None will be 
+                          returned (insitu=False). 
+                          (default is "": no table output)
 
         Example:
@@ -3082 +4219 @@
                 workscan = self.copy()
 
-            #if mask         is None: mask         = [True for i in xrange(workscan.nchan())]
             if mask         is None: mask         = []
             if order        is None: order        = 5
@@ -3093 +4229 @@
             if outlog       is None: outlog       = False
             if blfile       is None: blfile       = ''
-            if csvformat     is None: csvformat     = False
-
-            if csvformat:
-                scsvformat = "T"
-            else:
-                scsvformat = "F"
+            if csvformat    is None: csvformat    = False
+            if bltable      is None: bltable      = ''
+
+            scsvformat = 'T' if csvformat else 'F'
 
             #CURRENTLY, PLOT=true UNAVAILABLE UNTIL cubic spline fitting is implemented as a fitter method. 
-            workscan._chebyshev_baseline(mask, order, clipthresh, clipniter, 
+            workscan._chebyshev_baseline(mask, order,
+                                         clipthresh, clipniter, 
                                          getresidual, 
                                          pack_progress_params(showprogress, 
                                                               minnrow),
-                                         outlog, scsvformat+blfile)
+                                         outlog, scsvformat+blfile, 
+                                         bltable)
             workscan._add_history("chebyshev_baseline", varlist)
-            
-            if insitu:
-                self._assign(workscan)
+
+            if bltable == '':
+                if insitu:
+                    self._assign(workscan)
+                else:
+                    return workscan
             else:
-                return workscan
+                if not insitu:
+                    return None
             
         except RuntimeError, e:
@@ -3117 +4257 @@
 
     @asaplog_post_dec
-    def auto_chebyshev_baseline(self, insitu=None, mask=None, order=None, 
+    def auto_chebyshev_baseline(self, mask=None, order=None, insitu=None, 
                               clipthresh=None, clipniter=None,
                               edge=None, threshold=None, chan_avg_limit=None, 
                               getresidual=None, plot=None,
                               showprogress=None, minnrow=None, outlog=None,
-                              blfile=None, csvformat=None):
+                              blfile=None, csvformat=None, bltable=None):
         """\
         Return a scan which has been baselined (all rows) by Chebyshev polynomials. 
@@ -3129 +4269 @@
 
         Parameters:
+            mask:           an optional mask retreived from scantable
+            order:          the maximum order of Chebyshev polynomial (default is 5)
             insitu:         if False a new scantable is returned.
                             Otherwise, the scaling is done in-situ
                             The default is taken from .asaprc (False)
-            mask:           an optional mask retreived from scantable
-            order:          the maximum order of Chebyshev polynomial (default is 5)
             clipthresh:     Clipping threshold. (default is 3.0, unit: sigma)
             clipniter:      maximum number of iteration of 'clipthresh'-sigma 
@@ -3173 +4313 @@
                             (default is "": no file/logger output)
             csvformat:      if True blfile is csv-formatted, default is False.
+            bltable:        name of a baseline table where fitting results
+                            (coefficients, rms, etc.) are to be written.
+                            if given, fitting results will NOT be output to
+                            scantable (insitu=True) or None will be 
+                            returned (insitu=False). 
+                            (default is "": no table output)
 
         Example:
@@ -3191 +4337 @@
                 workscan = self.copy()
             
-            #if mask           is None: mask           = [True for i in xrange(workscan.nchan())]
             if mask           is None: mask           = []
             if order          is None: order          = 5
@@ -3206 +4351 @@
             if blfile         is None: blfile         = ''
             if csvformat      is None: csvformat      = False
-
-            if csvformat:
-                scsvformat = "T"
-            else:
-                scsvformat = "F"
+            if bltable        is None: bltable        = ''
+
+            scsvformat = 'T' if csvformat else 'F'
 
             #CURRENTLY, PLOT=true UNAVAILABLE UNTIL cubic spline fitting is implemented as a fitter method.
-            workscan._auto_chebyshev_baseline(mask, order, clipthresh, 
-                                              clipniter, 
+            workscan._auto_chebyshev_baseline(mask, order,
+                                              clipthresh, clipniter, 
                                               normalise_edge_param(edge), 
                                               threshold, 
@@ -3220 +4363 @@
                                               pack_progress_params(showprogress, 
                                                                    minnrow), 
-                                              outlog, scsvformat+blfile)
+                                              outlog, scsvformat+blfile,
+                                              bltable)
             workscan._add_history("auto_chebyshev_baseline", varlist)
-            
-            if insitu:
-                self._assign(workscan)
+
+            if bltable == '':
+                if insitu:
+                    self._assign(workscan)
+                else:
+                    return workscan
             else:
-                return workscan
+                if not insitu:
+                    return None
             
         except RuntimeError, e:
@@ -3232 +4380 @@
 
     @asaplog_post_dec
-    def poly_baseline(self, mask=None, order=None, insitu=None, plot=None, 
+    def poly_baseline(self, mask=None, order=None, insitu=None, 
+                      clipthresh=None, clipniter=None, plot=None, 
                       getresidual=None, showprogress=None, minnrow=None, 
-                      outlog=None, blfile=None, csvformat=None):
+                      outlog=None, blfile=None, csvformat=None,
+                      bltable=None):
         """\
         Return a scan which has been baselined (all rows) by a polynomial.
         Parameters:
+            mask:         an optional mask
+            order:        the order of the polynomial (default is 0)
             insitu:       if False a new scantable is returned.
                           Otherwise, the scaling is done in-situ
                           The default is taken from .asaprc (False)
-            mask:         an optional mask
-            order:        the order of the polynomial (default is 0)
+            clipthresh:   Clipping threshold. (default is 3.0, unit: sigma)
+            clipniter:    maximum number of iteration of 'clipthresh'-sigma 
+                          clipping (default is 0)
             plot:         plot the fit and the residual. In this each
                           indivual fit has to be approved, by typing 'y'
@@ -3258 +4411 @@
                           (default is "": no file/logger output)
             csvformat:    if True blfile is csv-formatted, default is False.
+            bltable:      name of a baseline table where fitting results
+                          (coefficients, rms, etc.) are to be written.
+                          if given, fitting results will NOT be output to
+                          scantable (insitu=True) or None will be 
+                          returned (insitu=False). 
+                          (default is "": no table output)
 
         Example:
@@ -3275 +4434 @@
                 workscan = self.copy()
 
-            #if mask         is None: mask         = [True for i in \
-            #                                           xrange(workscan.nchan())]
             if mask         is None: mask         = []
             if order        is None: order        = 0
+            if clipthresh   is None: clipthresh   = 3.0
+            if clipniter    is None: clipniter    = 0
             if plot         is None: plot         = False
             if getresidual  is None: getresidual  = True
@@ -3284 +4443 @@
             if minnrow      is None: minnrow      = 1000
             if outlog       is None: outlog       = False
-            if blfile       is None: blfile       = ""
+            if blfile       is None: blfile       = ''
             if csvformat    is None: csvformat    = False
-
-            if csvformat:
-                scsvformat = "T"
-            else:
-                scsvformat = "F"
+            if bltable      is None: bltable      = ''
+
+            scsvformat = 'T' if csvformat else 'F'
 
             if plot:
@@ -3354 +4511 @@
                     blf.close()
             else:
-                workscan._poly_baseline(mask, order, getresidual, 
+                workscan._poly_baseline(mask, order,
+                                        clipthresh, clipniter, #
+                                        getresidual, 
                                         pack_progress_params(showprogress, 
                                                              minnrow), 
-                                        outlog, scsvformat+blfile)
+                                        outlog, scsvformat+blfile,
+                                        bltable)  #
             
             workscan._add_history("poly_baseline", varlist)
@@ -3370 +4530 @@
 
     @asaplog_post_dec
-    def auto_poly_baseline(self, mask=None, order=None, edge=None, 
-                           threshold=None, chan_avg_limit=None,
-                           plot=None, insitu=None,
-                           getresidual=None, showprogress=None, 
-                           minnrow=None, outlog=None, blfile=None, csvformat=None):
+    def auto_poly_baseline(self, mask=None, order=None, insitu=None, 
+                           clipthresh=None, clipniter=None, 
+                           edge=None, threshold=None, chan_avg_limit=None,
+                           getresidual=None, plot=None, 
+                           showprogress=None, minnrow=None, outlog=None,
+                           blfile=None, csvformat=None, bltable=None):
         """\
         Return a scan which has been baselined (all rows) by a polynomial.
@@ -3381 +4542 @@
 
         Parameters:
+            mask:           an optional mask retreived from scantable
+            order:          the order of the polynomial (default is 0)
             insitu:         if False a new scantable is returned.
                             Otherwise, the scaling is done in-situ
                             The default is taken from .asaprc (False)
-            mask:           an optional mask retreived from scantable
-            order:          the order of the polynomial (default is 0)
+            clipthresh:     Clipping threshold. (default is 3.0, unit: sigma)
+            clipniter:      maximum number of iteration of 'clipthresh'-sigma 
+                            clipping (default is 0)
             edge:           an optional number of channel to drop at
                             the edge of spectrum. If only one value is
@@ -3421 +4585 @@
                             (default is "": no file/logger output)
             csvformat:      if True blfile is csv-formatted, default is False.
+            bltable:        name of a baseline table where fitting results
+                            (coefficients, rms, etc.) are to be written.
+                            if given, fitting results will NOT be output to
+                            scantable (insitu=True) or None will be 
+                            returned (insitu=False). 
+                            (default is "": no table output)
 
         Example:
@@ -3436 +4606 @@
                 workscan = self.copy()
 
-            #if mask           is None: mask           = [True for i in xrange(workscan.nchan())]
             if mask           is None: mask           = []
             if order          is None: order          = 0
+            if clipthresh     is None: clipthresh     = 3.0
+            if clipniter      is None: clipniter      = 0
             if edge           is None: edge           = (0, 0)
             if threshold      is None: threshold      = 3
@@ -3449 +4620 @@
             if blfile         is None: blfile         = ''
             if csvformat      is None: csvformat      = False
-
-            if csvformat:
-                scsvformat = "T"
-            else:
-                scsvformat = "F"
+            if bltable        is None: bltable        = ''
+
+            scsvformat = 'T' if csvformat else 'F'
 
             edge = normalise_edge_param(edge)
@@ -3523 +4692 @@
                 if outblfile: blf.close()
             else:
-                workscan._auto_poly_baseline(mask, order, edge, threshold, 
+                workscan._auto_poly_baseline(mask, order,
+                                             clipthresh, clipniter, 
+                                             edge, threshold, 
                                              chan_avg_limit, getresidual, 
                                              pack_progress_params(showprogress,
                                                                   minnrow), 
-                                             outlog, scsvformat+blfile)
-
+                                             outlog, scsvformat+blfile,
+                                             bltable)
             workscan._add_history("auto_poly_baseline", varlist)
-            
-            if insitu:
-                self._assign(workscan)
+
+            if bltable == '':
+                if insitu:
+                    self._assign(workscan)
+                else:
+                    return workscan
             else:
-                return workscan
+                if not insitu:
+                    return None
             
         except RuntimeError, e:
@@ -3644 +4819 @@
         self._math._setinsitu(insitu)
         varlist = vars()
-        s = scantable(self._math._unaryop(self, offset, "ADD", False))
+        s = scantable(self._math._unaryop(self, offset, "ADD", False, False))
         s._add_history("add", varlist)
         if insitu:
@@ -3667 +4842 @@
             tsys:        if True (default) then apply the operation to Tsys
                          as well as the data
-
         """
         if insitu is None: insitu = rcParams['insitu']
@@ -3678 +4852 @@
                                                          numpy.ndarray):
                 from asapmath import _array2dOp
-                s = _array2dOp( self, factor, "MUL", tsys, insitu )
+                s = _array2dOp( self, factor, "MUL", tsys, insitu, True )
             else:
                 s = scantable( self._math._arrayop( self, factor, 
-                                                    "MUL", tsys ) )
+                                                    "MUL", tsys, True ) )
         else:
-            s = scantable(self._math._unaryop(self, factor, "MUL", tsys))
+            s = scantable(self._math._unaryop(self, factor, "MUL", tsys, True ))
         s._add_history("scale", varlist)
         if insitu:
@@ -3730 +4904 @@
         self._add_history("set_sourcetype", varlist)
 
+
+    def set_sourcename(self, name):
+        varlist = vars()
+        self._setsourcename(name)
+        self._add_history("set_sourcename", varlist)
+
     @asaplog_post_dec
     @preserve_selection
@@ -3766 +4946 @@
         s = None
         if mode.lower() == "paired":
+            from asap._asap import srctype
             sel = self.get_selection()
-            sel.set_query("SRCTYPE==psoff")
+            #sel.set_query("SRCTYPE==psoff")
+            sel.set_types(srctype.psoff)
             self.set_selection(sel)
             offs = self.copy()
-            sel.set_query("SRCTYPE==pson")
+            #sel.set_query("SRCTYPE==pson")
+            sel.set_types(srctype.pson)
             self.set_selection(sel)
             ons = self.copy()
@@ -3887 +5070 @@
             if other == 0.0:
                 raise ZeroDivisionError("Dividing by zero is not recommended")
-            s = scantable(self._math._unaryop(self, other, mode, False))
+            s = scantable(self._math._unaryop(self, other, mode, False, True))
         elif isinstance(other, list) or isinstance(other, numpy.ndarray):
             if isinstance(other[0], list) \
                     or isinstance(other[0], numpy.ndarray):
                 from asapmath import _array2dOp
-                s = _array2dOp( self, other, mode, False )
+                s = _array2dOp(self, other, mode, False)
             else:
-                s = scantable( self._math._arrayop( self, other, 
-                                                    mode, False ) )
+                s = scantable(self._math._arrayop(self, other, mode, False, True))
         else:
             raise TypeError("Other input is not a scantable or float value")
@@ -3928 +5110 @@
             self.set_selection(basesel+sel)
             nans = numpy.isnan(self._getspectrum(0))
-        if numpy.any(nans):
-            bnans = [ bool(v) for v in nans]
-            self.flag(bnans)
+            if numpy.any(nans):
+                bnans = [ bool(v) for v in nans]
+                self.flag(bnans)
+        
+        self.set_selection(basesel)
 
     def get_row_selector(self, rowno):