Changeset 1446 for branches/alma/python


Ignore:
Timestamp:
11/12/08 17:04:01 (16 years ago)
Author:
TakTsutsumi
Message:

Merged recent updates (since 2007) from nrao-asap

Location:
branches/alma/python
Files:
5 edited

Legend:

Unmodified
Added
Removed

  • branches/alma/python/asapfitter.py

    r1389 r1446  
    586586            scan = self.data
    587587        rows = xrange(scan.nrow())
     588        # Save parameters of baseline fits as a class attribute.
     589        # NOTICE: This does not reflect changes in scantable!
     590        if len(rows) > 0: self.blpars=[]
    588591        from asap import asaplog
    589592        asaplog.push("Fitting:")
    … …  
    595598            self.data = None
    596599            self.fit()
    597             x = self.get_parameters()
     600            fpar = self.get_parameters()
    598601            if plot:
    599602                self.plot(residual=True)
    600603                x = raw_input("Accept fit ([y]/n): ")
    601604                if x.upper() == 'N':
     605                    self.blpars.append(None)
    602606                    continue
    603607            scan._setspectrum(self.fitter.getresidual(), r)
     608            self.blpars.append(fpar)
    604609        if plot:
    605610            self._p.unmap()

  • branches/alma/python/asaplotbase.py

    r1259 r1446  
    161161        y2 = range(l2)
    162162        m2 = range(l2)
    163         ymsk = y.raw_mask()
    164         ydat = y.raw_data()
     163        #ymsk = y.raw_mask()
     164        #ydat = y.raw_data()
     165        ymsk = y.mask
     166        ydat = y.data
    165167        for i in range(l2):
    166168            x2[i] = x[i/2]
    … …  
    707709            for sp in self.subplots:
    708710                ax = sp['axes']
    709                 s = ax.title.get_size()
    710                 tsize = s-(self.cols+self.rows)
    711                 ax.title.set_size(tsize)
     711                s = rcParams['axes.titlesize']
     712                tsize = s-(self.cols+self.rows-1)
     713                ax.title.set_size(max(tsize,9))
    712714                fp = FP(size=rcParams['axes.labelsize'])
    713715                setp(ax.get_xticklabels(), fontsize=xts)

  • branches/alma/python/asapmath.py

    r1389 r1446  
    4545    if kwargs.has_key('align'):
    4646        align = kwargs.get('align')
     47    compel = False
     48    if kwargs.has_key('compel'):
     49        compel = kwargs.get('compel')
    4750    varlist = vars()
    4851    if isinstance(args[0],list):
    … …  
    8689        del merged
    8790    else:
    88         s = scantable(stm._average(alignedlst, mask, weight.upper(), scanav))
     91        #s = scantable(stm._average(alignedlst, mask, weight.upper(), scanav))
     92        s = scantable(stm._new_average(alignedlst, compel, mask, weight.upper(), scanav))
    8993    s._add_history("average_time",varlist)
    9094    print_log()

  • branches/alma/python/asapplotter.py

    r1389 r1446  
    605605        if n > 1:
    606606            ganged = rcParams['plotter.ganged']
     607            ###Start Mod: 2008.09.22 kana ###
     608            if self._panelling == 'i':
     609                ganged = False
     610            ###End Mod#######################
    607611            if self._rows and self._cols:
    608612                n = min(n,self._rows*self._cols)
    … …  
    688692                if stackcount == nstack:
    689693                    allxlim.sort()
    690                     self._plotter.axes.set_xlim([allxlim[0],allxlim[-1]])
     694                    ###Start Mod: 2008.09.22 kana ###
     695                    #self._plotter.axes.set_xlim([allxlim[0],allxlim[-1]])
     696                    self._plotter.subplots[panelcount-1]['axes'].set_xlim([allxlim[0],allxlim[-1]])
     697                    ###End Mod#######################
    691698                    # clear
    692699                    allxlim =[]
    … …  
    741748        return userlabel or d[mode]
    742749
    743     def plotazel(self, scan=None):
     750    def plotazel(self, scan=None, outfile=None):
    744751        """
    745752        plot azimuth and elevation  versus time of a scantable
    … …  
    750757        from matplotlib.numerix import array, pi
    751758        self._data = scan
     759        self._outfile = outfile
    752760        dates = self._data.get_time()
    753761        t = PL.date2num(dates)
    754762        tz = timezone('UTC')
    755763        PL.cla()
    756         PL.ioff()
     764        #PL.ioff()
    757765        PL.clf()
    758766        tdel = max(t) - min(t)
    … …  
    772780            minloc = MinuteLocator(20)
    773781        PL.title(dstr)
    774         PL.plot_date(t,el,'b,', tz=tz)
    775         #ax.grid(True)
     782
     783        if tdel == 0.0:
     784            th = (t - PL.floor(t))*24.0
     785            PL.plot(th,el,'o',markersize=2, markerfacecolor='b', markeredgecolor='b')
     786        else:
     787            PL.plot_date(t,el,'o', markersize=2, markerfacecolor='b', markeredgecolor='b',tz=tz)
     788            #ax.grid(True)
     789            ax.xaxis.set_major_formatter(timefmt)
     790            ax.xaxis.set_major_locator(majloc)
     791            ax.xaxis.set_minor_locator(minloc)
    776792        ax.yaxis.grid(True)
    777793        yloc = MultipleLocator(30)
    778794        ax.set_ylim(0,90)
    779         ax.xaxis.set_major_formatter(timefmt)
    780         ax.xaxis.set_major_locator(majloc)
    781         ax.xaxis.set_minor_locator(minloc)
    782795        ax.yaxis.set_major_locator(yloc)
    783796        if tdel > 1.0:
    … …  
    785798        #    PL.setp(labels, fontsize=10, rotation=45)
    786799            PL.setp(labels, fontsize=10)
     800
    787801        # Az plot
    788802        az = array(self._data.get_azimuth())*180./pi
    … …  
    792806
    793807        ax = PL.subplot(2,1,2)
    794         PL.xlabel('Time (UT)')
     808        #PL.xlabel('Time (UT [hour])')
    795809        PL.ylabel('Az [deg.]')
    796         PL.plot_date(t,az,'b,', tz=tz)
     810        if tdel == 0.0:
     811            PL.plot(th,az,'o',markersize=2, markeredgecolor='b',markerfacecolor='b')
     812        else:
     813            PL.plot_date(t,az,'o', markersize=2,markeredgecolor='b',markerfacecolor='b',tz=tz)
     814            ax.xaxis.set_major_formatter(timefmt)
     815            ax.xaxis.set_major_locator(majloc)
     816            ax.xaxis.set_minor_locator(minloc)
     817        #ax.grid(True)
    797818        ax.set_ylim(0,360)
    798         #ax.grid(True)
    799819        ax.yaxis.grid(True)
    800820        #hfmt = DateFormatter('%H')
    801821        #hloc = HourLocator()
    802822        yloc = MultipleLocator(60)
    803         ax.xaxis.set_major_formatter(timefmt)
    804         ax.xaxis.set_major_locator(majloc)
    805         ax.xaxis.set_minor_locator(minloc)
    806823        ax.yaxis.set_major_locator(yloc)
    807824        if tdel > 1.0:
    808825            labels = ax.get_xticklabels()
    809826            PL.setp(labels, fontsize=10)
    810         PL.ion()
     827            PL.xlabel('Time (UT [day])')
     828        else:
     829            PL.xlabel('Time (UT [hour])')
     830
     831        #PL.ion()
    811832        PL.draw()
    812 
    813     def plotpointing(self, scan=None):
     833        if (self._outfile is not None):
     834           PL.savefig(self._outfile)
     835
     836    def plotpointing(self, scan=None, outfile=None):
    814837        """
    815838        plot telescope pointings
    … …  
    820843        from matplotlib.numerix import array, pi, zeros
    821844        self._data = scan
     845        self._outfile = outfile
    822846        dir = array(self._data.get_directionval()).transpose()
    823847        ra = dir[0]*180./pi
    824848        dec = dir[1]*180./pi
    825849        PL.cla()
    826         PL.ioff()
     850        #PL.ioff()
    827851        PL.clf()
    828852        ax = PL.axes([0.1,0.1,0.8,0.8])
    … …  
    835859        [xmin,xmax,ymin,ymax] = PL.axis()
    836860        PL.axis([xmax,xmin,ymin,ymax])
    837         PL.ion()
     861        #PL.ion()
    838862        PL.draw()
    839 
     863        if (self._outfile is not None):
     864           PL.savefig(self._outfile)
     865
     866    # plot total power data
     867    # plotting in time is not yet implemented..
     868    def plottp(self, scan=None, outfile=None):
     869        if self._plotter.is_dead:
     870            self._plotter = self._newplotter()
     871        self._plotter.hold()
     872        self._plotter.clear()
     873        from asap import scantable
     874        if not self._data and not scan:
     875            msg = "Input is not a scantable"
     876            if rcParams['verbose']:
     877                print msg
     878                return
     879            raise TypeError(msg)
     880        if isinstance(scan, scantable):
     881            if self._data is not None:
     882                if scan != self._data:
     883                    self._data = scan
     884                    # reset
     885                    self._reset()
     886            else:
     887                self._data = scan
     888                self._reset()
     889        # ranges become invalid when abcissa changes?
     890        #if self._abcunit and self._abcunit != self._data.get_unit():
     891        #    self._minmaxx = None
     892        #    self._minmaxy = None
     893        #    self._abcunit = self._data.get_unit()
     894        #    self._datamask = None
     895        self._plottp(self._data)
     896        if self._minmaxy is not None:
     897            self._plotter.set_limits(ylim=self._minmaxy)
     898        self._plotter.release()
     899        self._plotter.tidy()
     900        self._plotter.show(hardrefresh=False)
     901        print_log()
     902        return
     903
     904    def _plottp(self,scan):
     905        """
     906        private method for plotting total power data
     907        """
     908        from matplotlib.numerix import ma, array, arange, logical_not
     909        r=0
     910        nr = scan.nrow()
     911        a0,b0 = -1,-1
     912        allxlim = []
     913        allylim = []
     914        y=[]
     915        self._plotter.set_panels()
     916        self._plotter.palette(0)
     917        #title
     918        #xlab = self._abcissa and self._abcissa[panelcount] \
     919        #       or scan._getabcissalabel()
     920        #ylab = self._ordinate and self._ordinate[panelcount] \
     921        #       or scan._get_ordinate_label()
     922        xlab = self._abcissa or 'row number' #or Time
     923        ylab = self._ordinate or scan._get_ordinate_label()
     924        self._plotter.set_axes('xlabel',xlab)
     925        self._plotter.set_axes('ylabel',ylab)
     926        lbl = self._get_label(scan, r, 's', self._title)
     927        if isinstance(lbl, list) or isinstance(lbl, tuple):
     928        #    if 0 <= panelcount < len(lbl):
     929        #        lbl = lbl[panelcount]
     930        #    else:
     931                # get default label
     932             lbl = self._get_label(scan, r, self._panelling, None)
     933        self._plotter.set_axes('title',lbl)
     934        y=array(scan._get_column(scan._getspectrum,-1))
     935        m = array(scan._get_column(scan._getmask,-1))
     936        y = ma.masked_array(y,mask=logical_not(array(m,copy=False)))
     937        x = arange(len(y))
     938        # try to handle spectral data somewhat...
     939        l,m = y.shape
     940        if m > 1:
     941            y=y.mean(axis=1)
     942        plotit = self._plotter.plot
     943        llbl = self._get_label(scan, r, self._stacking, None)
     944        self._plotter.set_line(label=llbl)
     945        if len(x) > 0:
     946            plotit(x,y)
     947
     948
     949    # forwards to matplotlib.Figure.text
     950    def figtext(self, *args, **kwargs):
     951        """
     952        Add text to figure at location x,y (relative 0-1 coords).
     953        This method forwards *args and **kwargs to a Matplotlib method,
     954        matplotlib.Figure.text.
     955        See the method help for detailed information.
     956        """
     957        self._plotter.text(*args, **kwargs)
     958    # end matplotlib.Figure.text forwarding function
     959

  • branches/alma/python/scantable.py

    r1401 r1446  
    300300        self._setselection(selection)
    301301
     302    def get_row(self, row=0, insitu=None):
     303        """
     304        Select a row in the scantable.
     305        Return a scantable with single row.
     306        Parameters:
     307            row: row no of integration, default is 0.
     308            insitu: if False a new scantable is returned.
     309                    Otherwise, the scaling is done in-situ
     310                    The default is taken from .asaprc (False)
     311        """
     312        if insitu is None: insitu = rcParams['insitu']
     313        if not insitu:
     314            workscan = self.copy()
     315        else:
     316            workscan = self
     317        # Select a row
     318        sel=selector()
     319        sel.set_scans([workscan.getscan(row)])
     320        sel.set_cycles([workscan.getcycle(row)])
     321        sel.set_beams([workscan.getbeam(row)])
     322        sel.set_ifs([workscan.getif(row)])
     323        sel.set_polarisations([workscan.getpol(row)])
     324        sel.set_name(workscan._getsourcename(row))
     325        workscan.set_selection(sel)
     326        if not workscan.nrow() == 1:
     327            msg = "Cloud not identify single row. %d rows selected."%(workscan.nrow())
     328            raise RuntimeError(msg)
     329        del sel
     330        if insitu:
     331            self._assign(workscan)
     332        else:
     333            return workscan
     334
    302335    def stats(self, stat='stddev', mask=None):
    303336        """
    … …  
    528561        return self._get_column(self._getdirectionvec, row)
    529562
    530 
    531563    def set_unit(self, unit='channel'):
    532564        """
    … …  
    768800        return msk
    769801
    770     def get_restfreqs(self):
     802    def get_masklist(self, mask=None, row=0):
     803        """
     804        Compute and return a list of mask windows, [min, max].
     805        Parameters:
     806            mask:       channel mask, created with create_mask.
     807            row:        calcutate the masklist using the specified row
     808                        for unit conversions, default is row=0
     809                        only necessary if frequency varies over rows.
     810        Returns:
     811            [min, max], [min2, max2], ...
     812                Pairs of start/end points (inclusive)specifying
     813                the masked regions
     814        """
     815        if not (isinstance(mask,list) or isinstance(mask, tuple)):
     816            raise TypeError("The mask should be list or tuple.")
     817        if len(mask) < 2:
     818            raise TypeError("The mask elements should be > 1")
     819        if self.nchan() != len(mask):
     820            msg = "Number of channels in scantable != number of mask elements"
     821            raise TypeError(msg)
     822        data = self._getabcissa(row)
     823        u = self._getcoordinfo()[0]
     824        if rcParams['verbose']:
     825            if u == "": u = "channel"
     826            msg = "The current mask window unit is %s" % u
     827            i = self._check_ifs()
     828            if not i:
     829                msg += "\nThis mask is only valid for IF=%d" % (self.getif(i))
     830            asaplog.push(msg)
     831        masklist=[]
     832        ist, ien = None, None
     833        ist, ien=self.get_mask_indices(mask)
     834        if ist is not None and ien is not None:
     835            for i in xrange(len(ist)):
     836                range=[data[ist[i]],data[ien[i]]]
     837                range.sort()
     838                masklist.append([range[0],range[1]])
     839        return masklist
     840
     841    def get_mask_indices(self, mask=None):
     842        """
     843        Compute and Return lists of mask start indices and mask end indices.
     844         Parameters:
     845            mask:       channel mask, created with create_mask.
     846        Returns:
     847            List of mask start indices and that of mask end indices,
     848            i.e., [istart1,istart2,....], [iend1,iend2,....].
     849        """
     850        if not (isinstance(mask,list) or isinstance(mask, tuple)):
     851            raise TypeError("The mask should be list or tuple.")
     852        if len(mask) < 2:
     853            raise TypeError("The mask elements should be > 1")
     854        istart=[]
     855        iend=[]
     856        if mask[0]: istart.append(0)
     857        for i in range(len(mask)-1):
     858            if not mask[i] and mask[i+1]:
     859                istart.append(i+1)
     860            elif mask[i] and not mask[i+1]:
     861                iend.append(i)
     862        if mask[len(mask)-1]: iend.append(len(mask)-1)
     863        if len(istart) != len(iend):
     864            raise RuntimeError("Numbers of mask start != mask end.")
     865        for i in range(len(istart)):
     866            if istart[i] > iend[i]:
     867                raise RuntimeError("Mask start index > mask end index")
     868                break
     869        return istart,iend
     870
     871#    def get_restfreqs(self):
     872#        """
     873#        Get the restfrequency(s) stored in this scantable.
     874#        The return value(s) are always of unit 'Hz'
     875#        Parameters:
     876#            none
     877#        Returns:
     878#            a list of doubles
     879#        """
     880#        return list(self._getrestfreqs())
     881
     882    def get_restfreqs(self, ids=None):
    771883        """
    772884        Get the restfrequency(s) stored in this scantable.
    773885        The return value(s) are always of unit 'Hz'
    774886        Parameters:
    775             none
     887            ids: (optional) a list of MOLECULE_ID for that restfrequency(s) to
     888                 be retrieved
    776889        Returns:
    777             a list of doubles
    778         """
    779         return list(self._getrestfreqs())
    780 
     890            dictionary containing ids and a list of doubles for each id
     891        """
     892        if ids is None:
     893            rfreqs={}
     894            idlist = self.getmolnos()
     895            for i in idlist:
     896                rfreqs[i]=list(self._getrestfreqs(i))
     897            return rfreqs
     898        else:
     899            if type(ids)==list or type(ids)==tuple:
     900                rfreqs={}
     901                for i in ids:
     902                    rfreqs[i]=list(self._getrestfreqs(i))
     903                return rfreqs
     904            else:
     905                return list(self._getrestfreqs(ids))
     906            #return list(self._getrestfreqs(ids))
    781907
    782908    def set_restfreqs(self, freqs=None, unit='Hz'):
    783909        """
     910        ********NEED TO BE UPDATED begin************
    784911        Set or replace the restfrequency specified and
    785912        If the 'freqs' argument holds a scalar,
    … …  
    793920        E.g. 'freqs=[1e9, 2e9]'  would mean IF 0 gets restfreq 1e9 and
    794921        IF 1 gets restfreq 2e9.
     922        ********NEED TO BE UPDATED end************
    795923        You can also specify the frequencies via a linecatalog/
    796924
    … …  
    800928
    801929        Example:
    802             # set the given restfrequency for the whole table
     930            # set the given restfrequency for the all currently selected IFs
    803931            scan.set_restfreqs(freqs=1.4e9)
    804             # If thee number of IFs in the data is >= 2 IF0 gets the first
    805             # value IF1 the second...
    806             scan.set_restfreqs(freqs=[1.4e9, 1.67e9])
     932            # set multiple restfrequencies to all the selected data
     933            scan.set_restfreqs(freqs=[1.4e9, 1.41e9, 1.42e9])
     934            # If the number of IFs in the data is >= 2 the IF0 gets the first
     935            # value IF1 the second... NOTE that freqs needs to be
     936            # specified in list of list (e.g. [[],[],...] ).
     937            scan.set_restfreqs(freqs=[[1.4e9],[1.67e9]])
    807938            #set the given restfrequency for the whole table (by name)
    808939            scan.set_restfreqs(freqs="OH1667")
    … …  
    824955        # simple  value
    825956        if isinstance(freqs, int) or isinstance(freqs, float):
    826             self._setrestfreqs(freqs, "",unit)
     957            # TT mod
     958            #self._setrestfreqs(freqs, "",unit)
     959            self._setrestfreqs([freqs], [""],unit)
    827960        # list of values
    828961        elif isinstance(freqs, list) or isinstance(freqs, tuple):
    829962            # list values are scalars
    830963            if isinstance(freqs[-1], int) or isinstance(freqs[-1], float):
     964                self._setrestfreqs(freqs, [""],unit)
     965            # list values are tuples, (value, name)
     966            elif isinstance(freqs[-1], dict):
     967                #sel = selector()
     968                #savesel = self._getselection()
     969                #iflist = self.getifnos()
     970                #for i in xrange(len(freqs)):
     971                #    sel.set_ifs(iflist[i])
     972                #    self._setselection(sel)
     973                #    self._setrestfreqs(freqs[i], "",unit)
     974                #self._setselection(savesel)
     975                self._setrestfreqs(freqs["value"],
     976                                   freqs["name"], "MHz")
     977            elif isinstance(freqs[-1], list) or isinstance(freqs[-1], tuple):
    831978                sel = selector()
    832979                savesel = self._getselection()
    833980                iflist = self.getifnos()
    834                 for i in xrange(len(freqs)):
    835                     sel.set_ifs(iflist[i])
    836                     self._setselection(sel)
    837                     self._setrestfreqs(freqs[i], "",unit)
    838                 self._setselection(savesel)
    839             # list values are tuples, (value, name)
    840             elif isinstance(freqs[-1], dict):
    841                 sel = selector()
    842                 savesel = self._getselection()
    843                 iflist = self.getifnos()
     981                if len(freqs)>len(iflist):
     982                    raise ValueError("number of elements in list of list exeeds the current IF selections")
    844983                for i in xrange(len(freqs)):
    845984                    sel.set_ifs(iflist[i])
    … …  
    852991            sel = selector()
    853992            savesel = self._getselection()
    854             iflist = self.getifnos()
    855993            for i in xrange(freqs.nrow()):
    856994                sel.set_ifs(iflist[i])
    … …  
    12541392            f = fitter()
    12551393            f.set_scan(self, mask)
    1256             #f.set_function(poly=order)
    12571394            if uselin:
    12581395                f.set_function(lpoly=order)
    … …  
    12601397                f.set_function(poly=order)
    12611398            s = f.auto_fit(insitu, plot=plot)
     1399            # Save parameters of baseline fits as a class attribute.
     1400            # NOTICE: It does not reflect changes in scantable!
     1401            self.blpars = f.blpars
    12621402            s._add_history("poly_baseline", varlist)
    12631403            print_log()
    … …  
    13551495
    13561496        rows = range(workscan.nrow())
     1497        # Save parameters of baseline fits & masklists as a class attribute.
     1498        # NOTICE: It does not reflect changes in scantable!
     1499        if len(rows) > 0:
     1500            self.blpars=[]
     1501            self.masklists=[]
    13571502        asaplog.push("Processing:")
    13581503        for r in rows:
    … …  
    13711516            # setup line finder
    13721517            fl.find_lines(r, mask, curedge)
     1518            outmask=fl.get_mask()
    13731519            f.set_scan(workscan, fl.get_mask())
    13741520            f.x = workscan._getabcissa(r)
    … …  
    13761522            f.data = None
    13771523            f.fit()
    1378             x = f.get_parameters()
     1524           
     1525            # Show mask list
     1526            masklist=workscan.get_masklist(fl.get_mask(),row=r)
     1527            msg = "mask range: "+str(masklist)
     1528            asaplog.push(msg, False)
     1529
     1530            fpar = f.get_parameters()
    13791531            if plot:
    13801532                f.plot(residual=True)
    13811533                x = raw_input("Accept fit ( [y]/n ): ")
    13821534                if x.upper() == 'N':
     1535                    self.blpars.append(None)
     1536                    self.masklists.append(None)
    13831537                    continue
    13841538            workscan._setspectrum(f.fitter.getresidual(), r)
     1539            self.blpars.append(fpar)
     1540            self.masklists.append(masklist)
    13851541        if plot:
    13861542            f._p.unmap()
    … …  
    17731929        if unit is not None:
    17741930            self.set_fluxunit(unit)
    1775         self.set_freqframe(rcParams['scantable.freqframe'])
    1776 
     1931        #self.set_freqframe(rcParams['scantable.freqframe'])
     1932
Note: See TracChangeset for help on using the changeset viewer.