source: trunk/python/__init__.py@ 1604

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

Added scantable.parallactify

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