source: trunk/python/__init__.py @ 1586

Last change on this file since 1586 was 1586, checked in by Malte Marquarding, 15 years ago

Ticket #165: have removed the hard-coding of parallactifying the data. NOTE THIS breaks the Table structure as I have moved the PARANGLE column from the main table into the FOCUS table. We need to have a new release. Also one needs to explicitly tell the scantable via rc or member function to enable parallactifying

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 24.1 KB
Line 
1"""
2This is the ATNF Single Dish Analysis package.
3
4"""
5import os,sys,shutil, platform
6
7# Set up AIPSPATH and first time use of asap i.e. ~/.asap/*
8plf = None
9if sys.platform == "linux2":
10    if platform.architecture()[0] == '64bit':
11        plf = 'linux_64b'
12    else:
13        plf = 'linux_gnu'
14elif sys.platform == 'darwin':
15    plf = 'darwin'
16else:
17    # Shouldn't happen - default to linux
18    plf = 'linux'
19asapdata = __path__[-1]
20# Allow user defined data location
21if os.environ.has_key("ASAPDATA"):
22    if os.path.exists(os.environ["ASAPDATA"]):
23        asapdata = os.environ["ASAPDATA"]
24# use AIPSPATH if defined and "data" dir present
25if not os.environ.has_key("AIPSPATH") or \
26        not os.path.exists(os.environ["AIPSPATH"].split()[0]+"/data"):
27    os.environ["AIPSPATH"] = "%s %s somwhere" % ( asapdata, plf)
28# set up user space
29userdir = os.environ["HOME"]+"/.asap"
30if not os.path.exists(userdir):
31    print 'First time ASAP use. Setting up ~/.asap'
32    os.mkdir(userdir)
33    shutil.copyfile(asapdata+"/data/ipythonrc-asap", userdir+"/ipythonrc-asap")
34    shutil.copyfile(asapdata+"/data/ipy_user_conf.py",
35                    userdir+"/ipy_user_conf.py")
36    f = file(userdir+"/asapuserfuncs.py", "w")
37    f.close()
38    f = file(userdir+"/ipythonrc", "w")
39    f.close()
40else:
41    # upgrade to support later ipython versions
42    if not os.path.exists(userdir+"/ipy_user_conf.py"):
43        shutil.copyfile(asapdata+"/data/ipy_user_conf.py",
44                        userdir+"/ipy_user_conf.py")
45
46# remove from namespace
47del asapdata, userdir, shutil, platform
48
49def _validate_bool(b):
50    'Convert b to a boolean or raise'
51    bl = b.lower()
52    if bl in ('f', 'no', 'false', '0', 0): return False
53    elif bl in ('t', 'yes', 'true', '1', 1): return True
54    else:
55        raise ValueError('Could not convert "%s" to boolean' % b)
56
57def _validate_int(s):
58    'convert s to int or raise'
59    try: return int(s)
60    except ValueError:
61        raise ValueError('Could not convert "%s" to int' % s)
62
63def _asap_fname():
64    """
65    Return the path to the rc file
66
67    Search order:
68
69     * current working dir
70     * environ var ASAPRC
71     * HOME/.asaprc
72
73    """
74
75    fname = os.path.join( os.getcwd(), '.asaprc')
76    if os.path.exists(fname): return fname
77
78    if os.environ.has_key('ASAPRC'):
79        path =  os.environ['ASAPRC']
80        if os.path.exists(path):
81            fname = os.path.join(path, '.asaprc')
82            if os.path.exists(fname):
83                return fname
84
85    if os.environ.has_key('HOME'):
86        home =  os.environ['HOME']
87        fname = os.path.join(home, '.asaprc')
88        if os.path.exists(fname):
89            return fname
90    return None
91
92
93defaultParams = {
94    # general
95    'verbose'             : [True, _validate_bool],
96    'useplotter'          : [True, _validate_bool],
97    'insitu'              : [True, _validate_bool],
98
99    # plotting
100    'plotter.gui'         : [True, _validate_bool],
101    'plotter.stacking'    : ['p', str],
102    'plotter.panelling'   : ['s', str],
103    'plotter.colours'     : ['', str],
104    'plotter.linestyles'  : ['', str],
105    'plotter.decimate'    : [False, _validate_bool],
106    'plotter.ganged'      : [True, _validate_bool],
107    'plotter.histogram'  : [False, _validate_bool],
108    'plotter.papertype'  : ['A4', str],
109    'plotter.axesformatting' : ['mpl', str],
110
111    # scantable
112    'scantable.save'      : ['ASAP', str],
113    'scantable.autoaverage'      : [True, _validate_bool],
114    'scantable.freqframe' : ['LSRK', str],  #default frequency frame
115    'scantable.verbosesummary'   : [False, _validate_bool],
116    'scantable.storage'   : ['memory', str],
117    'scantable.history'   : [True, _validate_bool],
118    'scantable.reference'      : ['.*(e|w|_R)$', str],
119    'scantable.parallactify'   : [False, _validate_bool]
120    # fitter
121    }
122
123def list_rcparameters():
124
125    print """
126# general
127# print verbose output
128verbose                    : True
129
130# preload a default plotter
131useplotter                 : True
132
133# apply operations on the input scantable or return new one
134insitu                     : True
135
136# plotting
137
138# do we want a GUI or plot to a file
139plotter.gui                : True
140
141# default mode for colour stacking
142plotter.stacking           : Pol
143
144# default mode for panelling
145plotter.panelling          : scan
146
147# push panels together, to share axis labels
148plotter.ganged             : True
149
150# decimate the number of points plotted by a factor of
151# nchan/1024
152plotter.decimate           : False
153
154# default colours/linestyles
155plotter.colours            :
156plotter.linestyles         :
157
158# enable/disable histogram plotting
159plotter.histogram          : False
160
161# ps paper type
162plotter.papertype          : A4
163
164# The formatting style of the xaxis
165plotter.axesformatting    : 'mpl' (default) or 'asap' (for old versions of matplotlib)
166
167# scantable
168
169# default storage of scantable ('memory'/'disk')
170scantable.storage          : memory
171
172# write history of each call to scantable
173scantable.history          : True
174
175# default ouput format when saving
176scantable.save             : ASAP
177
178# auto averaging on read
179scantable.autoaverage      : True
180
181# default frequency frame to set when function
182# scantable.set_freqframe is called
183scantable.freqframe        : LSRK
184
185# Control the level of information printed by summary
186scantable.verbosesummary   : False
187
188# Control the identification of reference (off) scans
189# This is has to be a regular expression
190scantable.reference         : .*(e|w|_R)$
191
192# Fitter
193"""
194
195def rc_params():
196    'Return the default params updated from the values in the rc file'
197
198    fname = _asap_fname()
199
200    if fname is None or not os.path.exists(fname):
201        message = 'could not find rc file; returning defaults'
202        ret =  dict([ (key, tup[0]) for key, tup in defaultParams.items()])
203        #print message
204        return ret
205
206    cnt = 0
207    for line in file(fname):
208        cnt +=1
209        line = line.strip()
210        if not len(line): continue
211        if line.startswith('#'): continue
212        tup = line.split(':',1)
213        if len(tup) !=2:
214            print ('Illegal line #%d\n\t%s\n\tin file "%s"' % (cnt, line, fname))
215            continue
216
217        key, val = tup
218        key = key.strip()
219        if not defaultParams.has_key(key):
220            print ('Bad key "%s" on line %d in %s' % (key, cnt, fname))
221            continue
222
223        default, converter =  defaultParams[key]
224
225        ind = val.find('#')
226        if ind>=0: val = val[:ind]   # ignore trailing comments
227        val = val.strip()
228        try: cval = converter(val)   # try to convert to proper type or raise
229        except ValueError, msg:
230            print ('Bad val "%s" on line #%d\n\t"%s"\n\tin file "%s"\n\t%s' % (val, cnt, line, fname, msg))
231            continue
232        else:
233            # Alles Klar, update dict
234            defaultParams[key][0] = cval
235
236    # strip the conveter funcs and return
237    ret =  dict([ (key, tup[0]) for key, tup in defaultParams.items()])
238    print ('loaded rc file %s'%fname)
239
240    return ret
241
242
243# this is the instance used by the asap classes
244rcParams = rc_params()
245
246rcParamsDefault = dict(rcParams.items()) # a copy
247
248def rc(group, **kwargs):
249    """
250    Set the current rc params.  Group is the grouping for the rc, eg
251    for scantable.save the group is 'scantable', for plotter.stacking, the
252    group is 'plotter', and so on.  kwargs is a list of attribute
253    name/value pairs, eg
254
255      rc('scantable', save='SDFITS')
256
257    sets the current rc params and is equivalent to
258
259      rcParams['scantable.save'] = 'SDFITS'
260
261    Use rcdefaults to restore the default rc params after changes.
262    """
263
264    aliases = {}
265
266    for k,v in kwargs.items():
267        name = aliases.get(k) or k
268        if len(group):
269            key = '%s.%s' % (group, name)
270        else:
271            key = name
272        if not rcParams.has_key(key):
273            raise KeyError('Unrecognized key "%s" for group "%s" and name "%s"' % (key, group, name))
274
275        rcParams[key] = v
276
277
278def rcdefaults():
279    """
280    Restore the default rc params - the ones that were created at
281    asap load time
282    """
283    rcParams.update(rcParamsDefault)
284
285def _n_bools(n, val):
286    return [ val for i in xrange(n) ]
287
288def _is_sequence_or_number(param, ptype=int):
289    if isinstance(param,tuple) or isinstance(param,list):
290        if len(param) == 0: return True # empty list
291        out = True
292        for p in param:
293            out &= isinstance(p,ptype)
294        return out
295    elif isinstance(param, ptype):
296        return True
297    return False
298
299def _to_list(param, ptype=int):
300    if isinstance(param, ptype):
301        if ptype is str: return param.split()
302        else: return [param]
303    if _is_sequence_or_number(param, ptype):
304        return param
305    return None
306
307def unique(x):
308    """
309    Return the unique values in a list
310    Parameters:
311        x:      the list to reduce
312    Examples:
313        x = [1,2,3,3,4]
314        print unique(x)
315        [1,2,3,4]
316    """
317    return dict([ (val, 1) for val in x]).keys()
318
319def list_files(path=".",suffix="rpf"):
320    """
321    Return a list files readable by asap, such as rpf, sdfits, mbf, asap
322    Parameters:
323        path:     The directory to list (default '.')
324        suffix:   The file extension (default rpf)
325    Example:
326        files = list_files("data/","sdfits")
327        print files
328        ['data/2001-09-01_0332_P363.sdfits',
329        'data/2003-04-04_131152_t0002.sdfits',
330        'data/Sgr_86p262_best_SPC.sdfits']
331    """
332    if not os.path.isdir(path):
333        return None
334    valid = "rpf rpf.1 rpf.2 sdf sdfits mbf asap".split()
335    if not suffix in valid:
336        return None
337    files = [os.path.expanduser(os.path.expandvars(path+"/"+f)) for f in os.listdir(path)]
338    return filter(lambda x: x.endswith(suffix),files)
339
340# workaround for ipython, which redirects this if banner=0 in ipythonrc
341sys.stdout = sys.__stdout__
342sys.stderr = sys.__stderr__
343
344# Logging
345from asap._asap import Log as _asaplog
346global asaplog
347asaplog=_asaplog()
348if rcParams['verbose']:
349    asaplog.enable()
350else:
351    asaplog.disable()
352
353def print_log():
354    log = asaplog.pop()
355    if len(log) and rcParams['verbose']: print log
356    return
357
358def mask_and(a, b):
359    assert(len(a)==len(b))
360    return [ a[i] & b[i] for i in xrange(len(a)) ]
361
362def mask_or(a, b):
363    assert(len(a)==len(b))
364    return [ a[i] | b[i] for i in xrange(len(a)) ]
365
366def mask_not(a):
367    return [ not i for i in a ]
368
369from asapfitter import fitter
370from asapreader import reader
371from selector import selector
372
373from asapmath import *
374from scantable import scantable
375from asaplinefind import linefinder
376from linecatalog import linecatalog
377
378if rcParams['useplotter']:
379    try:
380        from asapplotter import asapplotter
381        gui = os.environ.has_key('DISPLAY') and rcParams['plotter.gui']
382        if gui:
383            import matplotlib
384            matplotlib.use("TkAgg")
385        import pylab
386        xyplotter = pylab
387        plotter = asapplotter(gui)
388        del gui
389    except ImportError:
390        print "Matplotlib not installed. No plotting available"
391
392__date__ = '$Date: 2009-06-30 05:46:14 +0000 (Tue, 30 Jun 2009) $'.split()[1]
393__version__  = '$Revision: 1586 $'
394
395def is_ipython():
396    return '__IP' in dir(sys.modules["__main__"])
397if is_ipython():
398    def version(): print  "ASAP %s(%s)"% (__version__, __date__)
399   
400    def list_scans(t = scantable):       
401        import inspect
402        print "The user created scantables are: ",
403        globs=inspect.currentframe().f_back.f_locals.copy()
404        out = [ k for k,v in globs.iteritems() \
405                     if isinstance(v, scantable) and not k.startswith("_") ]
406        print out
407        return out
408
409    def commands():
410        x = """
411    [The scan container]
412        scantable           - a container for integrations/scans
413                              (can open asap/rpfits/sdfits and ms files)
414            copy            - returns a copy of a scan
415            get_scan        - gets a specific scan out of a scantable
416                              (by name or number)
417            drop_scan       - drops a specific scan out of a scantable
418                              (by number)
419            set_selection   - set a new subselection of the data
420            get_selection   - get the current selection object
421            summary         - print info about the scantable contents
422            stats           - get specified statistic of the spectra in
423                              the scantable
424            stddev          - get the standard deviation of the spectra
425                              in the scantable
426            get_tsys        - get the TSys
427            get_time        - get the timestamps of the integrations
428            get_inttime     - get the integration time
429            get_sourcename  - get the source names of the scans
430            get_azimuth     - get the azimuth of the scans
431            get_elevation   - get the elevation of the scans
432            get_parangle    - get the parallactic angle of the scans
433            get_unit        - get the current unit
434            set_unit        - set the abcissa unit to be used from this
435                              point on
436            get_abcissa     - get the abcissa values and name for a given
437                              row (time)
438            get_column_names - get the names of the columns in the scantable
439                               for use with selector.set_query
440            set_freqframe   - set the frame info for the Spectral Axis
441                              (e.g. 'LSRK')
442            set_doppler     - set the doppler to be used from this point on
443            set_dirframe    - set the frame for the direction on the sky
444            set_instrument  - set the instrument name
445            set_feedtype    - set the feed type
446            get_fluxunit    - get the brightness flux unit
447            set_fluxunit    - set the brightness flux unit
448            set_sourcetype  - set the type of the source - source or reference
449            create_mask     - return an mask in the current unit
450                              for the given region. The specified regions
451                              are NOT masked
452            get_restfreqs   - get the current list of rest frequencies
453            set_restfreqs   - set a list of rest frequencies
454            shift_refpix    - shift the reference pixel of the IFs
455            set_spectrum    - overwrite the spectrum for a given row
456            get_spectrum    - retrieve the spectrum for a given
457            get_mask        - retrieve the mask for a given
458            flag            - flag selected channels in the data
459            lag_flag        - flag specified frequency in the data
460            save            - save the scantable to disk as either 'ASAP',
461                              'SDFITS' or 'ASCII'
462            nbeam,nif,nchan,npol - the number of beams/IFs/Pols/Chans
463            nscan           - the number of scans in the scantable
464            nrow            - the number of spectra in the scantable
465            history         - print the history of the scantable
466            get_fit         - get a fit which has been stored witnh the data
467            average_time    - return the (weighted) time average of a scan
468                              or a list of scans
469            average_pol     - average the polarisations together.
470            average_beam    - average the beams together.
471            convert_pol     - convert to a different polarisation type
472            auto_quotient   - return the on/off quotient with
473                              automatic detection of the on/off scans (closest
474                              in time off is selected)
475            mx_quotient     - Form a quotient using MX data (off beams)
476            scale, *, /     - return a scan scaled by a given factor
477            add, +          - return a scan with given value added
478            sub, -          - return a scan with given value subtracted
479            bin             - return a scan with binned channels
480            resample        - return a scan with resampled channels
481            smooth          - return the spectrally smoothed scan
482            poly_baseline   - fit a polynomial baseline to all Beams/IFs/Pols
483            auto_poly_baseline - automatically fit a polynomial baseline
484            recalc_azel     - recalculate azimuth and elevation based on
485                              the pointing
486            gain_el         - apply gain-elevation correction
487            opacity         - apply opacity correction
488            convert_flux    - convert to and from Jy and Kelvin brightness
489                              units
490            freq_align      - align spectra in frequency frame
491            invert_phase    - Invert the phase of the cross-correlation
492            swap_linears    - Swap XX and YY (or RR LL)
493            rotate_xyphase  - rotate XY phase of cross correlation
494            rotate_linpolphase - rotate the phase of the complex
495                                 polarization O=Q+iU correlation
496            freq_switch     - perform frequency switching on the data
497            stats           - Determine the specified statistic, e.g. 'min'
498                              'max', 'rms' etc.
499            stddev          - Determine the standard deviation of the current
500                              beam/if/pol
501     [Selection]
502         selector              - a selection object to set a subset of a scantable
503            set_scans          - set (a list of) scans by index
504            set_cycles         - set (a list of) cycles by index
505            set_beams          - set (a list of) beamss by index
506            set_ifs            - set (a list of) ifs by index
507            set_polarisations  - set (a list of) polarisations by name
508                                 or by index
509            set_names          - set a selection by name (wildcards allowed)
510            set_tsys           - set a selection by tsys thresholds
511            set_query          - set a selection by SQL-like query, e.g. BEAMNO==1
512            ( also  get_ functions for all these )
513            reset              - unset all selections
514            +                  - merge two selections
515
516     [Math] Mainly functions which operate on more than one scantable
517
518            average_time    - return the (weighted) time average
519                              of a list of scans
520            quotient        - return the on/off quotient
521            simple_math     - simple mathematical operations on two scantables,
522                              'add', 'sub', 'mul', 'div'
523            quotient        - build quotient of the given on and off scans
524                              (matched pairs and 1 off - n on are valid)
525            merge           - merge a list of scantables
526
527     [Line Catalog]
528        linecatalog              - a linecatalog wrapper, taking an ASCII or
529                                   internal format table
530            summary              - print a summary of the current selection
531            set_name             - select a subset by name pattern, e.g. '*OH*'
532            set_strength_limits  - select a subset by line strength limits
533            set_frequency_limits - select a subset by frequency limits
534            reset                - unset all selections
535            save                 - save the current subset to a table (internal
536                                   format)
537            get_row              - get the name and frequency from a specific
538                                   row in the table
539     [Fitting]
540        fitter
541            auto_fit        - return a scan where the function is
542                              applied to all Beams/IFs/Pols.
543            commit          - return a new scan where the fits have been
544                              commited.
545            fit             - execute the actual fitting process
546            store_fit       - store the fit parameters in the data (scantable)
547            get_chi2        - get the Chi^2
548            set_scan        - set the scantable to be fit
549            set_function    - set the fitting function
550            set_parameters  - set the parameters for the function(s), and
551                              set if they should be held fixed during fitting
552            set_gauss_parameters - same as above but specialised for individual
553                                   gaussian components
554            get_parameters  - get the fitted parameters
555            plot            - plot the resulting fit and/or components and
556                              residual
557    [Plotter]
558        asapplotter         - a plotter for asap, default plotter is
559                              called 'plotter'
560            plot            - plot a scantable
561            plot_lines      - plot a linecatalog overlay
562            save            - save the plot to a file ('png' ,'ps' or 'eps')
563            set_mode        - set the state of the plotter, i.e.
564                              what is to be plotted 'colour stacked'
565                              and what 'panelled'
566            set_selection   - only plot a selected part of the data
567            set_range       - set a 'zoom' window [xmin,xmax,ymin,ymax]
568            set_legend      - specify user labels for the legend indeces
569            set_title       - specify user labels for the panel indeces
570            set_abcissa     - specify a user label for the abcissa
571            set_ordinate    - specify a user label for the ordinate
572            set_layout      - specify the multi-panel layout (rows,cols)
573            set_colors      - specify a set of colours to use
574            set_linestyles  - specify a set of linestyles to use if only
575                              using one color
576            set_font        - set general font properties, e.g. 'family'
577            set_histogram   - plot in historam style
578            set_mask        - set a plotting mask for a specific polarization
579            text            - draw text annotations either in data or relative
580                              coordinates
581            arrow           - draw arrow annotations either in data or relative
582                              coordinates
583            axhline,axvline - draw horizontal/vertical lines
584            axhspan,axvspan - draw horizontal/vertical regions
585            annotate        - draw an arrow with label
586            create_mask     - create a scnatble mask interactively
587
588        xyplotter           - matplotlib/pylab plotting functions
589
590    [Reading files]
591        reader              - access rpfits/sdfits files
592            open            - attach reader to a file
593            close           - detach reader from file
594            read            - read in integrations
595            summary         - list info about all integrations
596
597    [General]
598        commands            - this command
599        print               - print details about a variable
600        list_scans          - list all scantables created by the user
601        list_files          - list all files readable by asap (default rpf)
602        del                 - delete the given variable from memory
603        range               - create a list of values, e.g.
604                              range(3) = [0,1,2], range(2,5) = [2,3,4]
605        help                - print help for one of the listed functions
606        execfile            - execute an asap script, e.g. execfile('myscript')
607        list_rcparameters   - print out a list of possible values to be
608                              put into $HOME/.asaprc
609        rc                  - set rc parameters from within asap
610        mask_and,mask_or,
611        mask_not            - boolean operations on masks created with
612                              scantable.create_mask
613
614    Note:
615        How to use this with help:
616                                         # function 'summary'
617        [xxx] is just a category
618        Every 'sub-level' in this list should be replaces by a '.' Period when
619        using help
620        Example:
621            ASAP> help scantable # to get info on ths scantable
622            ASAP> help scantable.summary # to get help on the scantable's
623            ASAP> help average_time
624
625            """
626        if rcParams['verbose']:
627            try:
628                from IPython.genutils import page as pager
629            except ImportError:
630                from pydoc import pager
631            pager(x)
632        else:
633            print x
634        return
635
636def welcome():
637    return """Welcome to ASAP v%s (%s) - the ATNF Spectral Analysis Package
638
639Please report any bugs via:
640http://svn.atnf.csiro.au/trac/asap/simpleticket
641
642[IMPORTANT: ASAP is 0-based]
643Type commands() to get a list of all available ASAP commands.""" % (__version__, __date__)
Note: See TracBrowser for help on using the repository browser.