source: trunk/python/__init__.py@ 1624

Last change on this file since 1624 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
Line 
1"""
2This is the ATNF Single Dish Analysis package.
3
4"""
5import os,sys,shutil, platform
6import functools
7
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]
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"]
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"):
28 os.environ["AIPSPATH"] = "%s %s somwhere" % ( asapdata, plf)
29# set up user space
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")
35 shutil.copyfile(asapdata+"/data/ipy_user_conf.py",
36 userdir+"/ipy_user_conf.py")
37 f = file(userdir+"/asapuserfuncs.py", "w")
38 f.close()
39 f = file(userdir+"/ipythonrc", "w")
40 f.close()
41else:
42 # upgrade to support later ipython versions
43 if not os.path.exists(userdir+"/ipy_user_conf.py"):
44 shutil.copyfile(asapdata+"/data/ipy_user_conf.py",
45 userdir+"/ipy_user_conf.py")
46
47# remove from namespace
48del asapdata, userdir, shutil, platform
49
50def _validate_bool(b):
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
58def _validate_int(s):
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
64def _asap_fname():
65 """
66 Return the path to the rc file
67
68 Search order:
69
70 * current working dir
71 * environ var ASAPRC
72 * HOME/.asaprc
73
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
93
94defaultParams = {
95 # general
96 'verbose' : [True, _validate_bool],
97 'useplotter' : [True, _validate_bool],
98 'insitu' : [True, _validate_bool],
99
100 # plotting
101 'plotter.gui' : [True, _validate_bool],
102 'plotter.stacking' : ['p', str],
103 'plotter.panelling' : ['s', str],
104 'plotter.colours' : ['', str],
105 'plotter.linestyles' : ['', str],
106 'plotter.decimate' : [False, _validate_bool],
107 'plotter.ganged' : [True, _validate_bool],
108 'plotter.histogram' : [False, _validate_bool],
109 'plotter.papertype' : ['A4', str],
110 'plotter.axesformatting' : ['mpl', str],
111
112 # scantable
113 'scantable.save' : ['ASAP', str],
114 'scantable.autoaverage' : [True, _validate_bool],
115 'scantable.freqframe' : ['LSRK', str], #default frequency frame
116 'scantable.verbosesummary' : [False, _validate_bool],
117 'scantable.storage' : ['memory', str],
118 'scantable.history' : [True, _validate_bool],
119 'scantable.reference' : ['.*(e|w|_R)$', str],
120 'scantable.parallactify' : [False, _validate_bool]
121 # fitter
122 }
123
124def list_rcparameters():
125
126 print """
127# general
128# print verbose output
129verbose : True
130
131# preload a default plotter
132useplotter : True
133
134# apply operations on the input scantable or return new one
135insitu : True
136
137# plotting
138
139# do we want a GUI or plot to a file
140plotter.gui : True
141
142# default mode for colour stacking
143plotter.stacking : Pol
144
145# default mode for panelling
146plotter.panelling : scan
147
148# push panels together, to share axis labels
149plotter.ganged : True
150
151# decimate the number of points plotted by a factor of
152# nchan/1024
153plotter.decimate : False
154
155# default colours/linestyles
156plotter.colours :
157plotter.linestyles :
158
159# enable/disable histogram plotting
160plotter.histogram : False
161
162# ps paper type
163plotter.papertype : A4
164
165# The formatting style of the xaxis
166plotter.axesformatting : 'mpl' (default) or 'asap' (for old versions of matplotlib)
167
168# scantable
169
170# default storage of scantable ('memory'/'disk')
171scantable.storage : memory
172
173# write history of each call to scantable
174scantable.history : True
175
176# default ouput format when saving
177scantable.save : ASAP
178
179# auto averaging on read
180scantable.autoaverage : True
181
182# default frequency frame to set when function
183# scantable.set_freqframe is called
184scantable.freqframe : LSRK
185
186# Control the level of information printed by summary
187scantable.verbosesummary : False
188
189# Control the identification of reference (off) scans
190# This is has to be a regular expression
191scantable.reference : .*(e|w|_R)$
192
193# Indicate whether the data was parallactified (total phase offest == 0.0)
194scantable.parallactify : False
195
196# Fitter
197"""
198
199def rc_params():
200 'Return the default params updated from the values in the rc file'
201
202 fname = _asap_fname()
203
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
209
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
220
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
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
233 except ValueError, msg:
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()])
242 print ('loaded rc file %s'%fname)
243
244 return ret
245
246
247# this is the instance used by the asap classes
248rcParams = rc_params()
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
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
257 name/value pairs, eg
258
259 rc('scantable', save='SDFITS')
260
261 sets the current rc params and is equivalent to
262
263 rcParams['scantable.save'] = 'SDFITS'
264
265 Use rcdefaults to restore the default rc params after changes.
266 """
267
268 aliases = {}
269
270 for k,v in kwargs.items():
271 name = aliases.get(k) or k
272 if len(group):
273 key = '%s.%s' % (group, name)
274 else:
275 key = name
276 if not rcParams.has_key(key):
277 raise KeyError('Unrecognized key "%s" for group "%s" and name "%s"' % (key, group, name))
278
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
289def _n_bools(n, val):
290 return [ val for i in xrange(n) ]
291
292def _is_sequence_or_number(param, ptype=int):
293 if isinstance(param,tuple) or isinstance(param,list):
294 if len(param) == 0: return True # empty list
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
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
310
311def unique(x):
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 """
321 return dict([ (val, 1) for val in x]).keys()
322
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
338 valid = "rpf rpf.1 rpf.2 sdf sdfits mbf asap".split()
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
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
351asaplog=_asaplog()
352if rcParams['verbose']:
353 asaplog.enable()
354else:
355 asaplog.disable()
356
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
366def print_log():
367 log = asaplog.pop().strip()
368 if len(log) and rcParams['verbose']: print log
369 return
370
371def mask_and(a, b):
372 assert(len(a)==len(b))
373 return [ a[i] & b[i] for i in xrange(len(a)) ]
374
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
382from asapfitter import fitter
383from asapreader import reader
384from selector import selector
385
386from asapmath import *
387from scantable import scantable
388from asaplinefind import linefinder
389from linecatalog import linecatalog
390
391if rcParams['useplotter']:
392 try:
393 from asapplotter import asapplotter
394 gui = os.environ.has_key('DISPLAY') and rcParams['plotter.gui']
395 if gui:
396 import matplotlib
397 matplotlib.use("TkAgg")
398 import pylab
399 xyplotter = pylab
400 plotter = asapplotter(gui)
401 del gui
402 except ImportError:
403 print "Matplotlib not installed. No plotting available"
404
405__date__ = '$Date: 2009-07-06 01:13:06 +0000 (Mon, 06 Jul 2009) $'.split()[1]
406__version__ = '$Revision: 1595 $'
407
408def is_ipython():
409 return '__IP' in dir(sys.modules["__main__"])
410
411if is_ipython():
412 def version(): print "ASAP %s(%s)"% (__version__, __date__)
413
414 def list_scans(t = scantable):
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
422
423 def commands():
424 x = """
425 [The scan container]
426 scantable - a container for integrations/scans
427 (can open asap/rpfits/sdfits and ms files)
428 copy - returns a copy of a scan
429 get_scan - gets a specific scan out of a scantable
430 (by name or number)
431 drop_scan - drops a specific scan out of a scantable
432 (by number)
433 set_selection - set a new subselection of the data
434 get_selection - get the current selection object
435 summary - print info about the scantable contents
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
440 get_tsys - get the TSys
441 get_time - get the timestamps of the integrations
442 get_inttime - get the integration time
443 get_sourcename - get the source names of the scans
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
447 get_unit - get the current unit
448 set_unit - set the abcissa unit to be used from this
449 point on
450 get_abcissa - get the abcissa values and name for a given
451 row (time)
452 get_column_names - get the names of the columns in the scantable
453 for use with selector.set_query
454 set_freqframe - set the frame info for the Spectral Axis
455 (e.g. 'LSRK')
456 set_doppler - set the doppler to be used from this point on
457 set_dirframe - set the frame for the direction on the sky
458 set_instrument - set the instrument name
459 set_feedtype - set the feed type
460 get_fluxunit - get the brightness flux unit
461 set_fluxunit - set the brightness flux unit
462 set_sourcetype - set the type of the source - source or reference
463 create_mask - return an mask in the current unit
464 for the given region. The specified regions
465 are NOT masked
466 get_restfreqs - get the current list of rest frequencies
467 set_restfreqs - set a list of rest frequencies
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
472 flag - flag selected channels in the data
473 lag_flag - flag specified frequency in the data
474 save - save the scantable to disk as either 'ASAP',
475 'SDFITS' or 'ASCII'
476 nbeam,nif,nchan,npol - the number of beams/IFs/Pols/Chans
477 nscan - the number of scans in the scantable
478 nrow - the number of spectra in the scantable
479 history - print the history of the scantable
480 get_fit - get a fit which has been stored witnh the data
481 average_time - return the (weighted) time average of a scan
482 or a list of scans
483 average_pol - average the polarisations together.
484 average_beam - average the beams together.
485 convert_pol - convert to a different polarisation type
486 auto_quotient - return the on/off quotient with
487 automatic detection of the on/off scans (closest
488 in time off is selected)
489 mx_quotient - Form a quotient using MX data (off beams)
490 scale, *, / - return a scan scaled by a given factor
491 add, + - return a scan with given value added
492 sub, - - return a scan with given value subtracted
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
497 auto_poly_baseline - automatically fit a polynomial baseline
498 recalc_azel - recalculate azimuth and elevation based on
499 the pointing
500 gain_el - apply gain-elevation correction
501 opacity - apply opacity correction
502 convert_flux - convert to and from Jy and Kelvin brightness
503 units
504 freq_align - align spectra in frequency frame
505 invert_phase - Invert the phase of the cross-correlation
506 swap_linears - Swap XX and YY (or RR LL)
507 rotate_xyphase - rotate XY phase of cross correlation
508 rotate_linpolphase - rotate the phase of the complex
509 polarization O=Q+iU correlation
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
515 get_row_selector - get the selection object for a specified row
516 number
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
527 set_query - set a selection by SQL-like query, e.g. BEAMNO==1
528 ( also get_ functions for all these )
529 reset - unset all selections
530 + - merge two selections
531
532 [Math] Mainly functions which operate on more than one scantable
533
534 average_time - return the (weighted) time average
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'
539 quotient - build quotient of the given on and off scans
540 (matched pairs and 1 off - n on are valid)
541 merge - merge a list of scantables
542
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
555 [Fitting]
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
562 store_fit - store the fit parameters in the data (scantable)
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
568 set_gauss_parameters - same as above but specialised for individual
569 gaussian components
570 get_parameters - get the fitted parameters
571 plot - plot the resulting fit and/or components and
572 residual
573 [Plotter]
574 asapplotter - a plotter for asap, default plotter is
575 called 'plotter'
576 plot - plot a scantable
577 plot_lines - plot a linecatalog overlay
578 save - save the plot to a file ('png' ,'ps' or 'eps')
579 set_mode - set the state of the plotter, i.e.
580 what is to be plotted 'colour stacked'
581 and what 'panelled'
582 set_selection - only plot a selected part of the data
583 set_range - set a 'zoom' window [xmin,xmax,ymin,ymax]
584 set_legend - specify user labels for the legend indeces
585 set_title - specify user labels for the panel indeces
586 set_abcissa - specify a user label for the abcissa
587 set_ordinate - specify a user label for the ordinate
588 set_layout - specify the multi-panel layout (rows,cols)
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
592 set_font - set general font properties, e.g. 'family'
593 set_histogram - plot in historam style
594 set_mask - set a plotting mask for a specific polarization
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
601 annotate - draw an arrow with label
602 create_mask - create a scnatble mask interactively
603
604 xyplotter - matplotlib/pylab plotting functions
605
606 [Reading files]
607 reader - access rpfits/sdfits files
608 open - attach reader to a file
609 close - detach reader from file
610 read - read in integrations
611 summary - list info about all integrations
612
613 [General]
614 commands - this command
615 print - print details about a variable
616 list_scans - list all scantables created by the user
617 list_files - list all files readable by asap (default rpf)
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')
623 list_rcparameters - print out a list of possible values to be
624 put into $HOME/.asaprc
625 rc - set rc parameters from within asap
626 mask_and,mask_or,
627 mask_not - boolean operations on masks created with
628 scantable.create_mask
629
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
635 using help
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
641 """
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
650 return
651
652def welcome():
653 return """Welcome to ASAP v%s (%s) - the ATNF Spectral Analysis Package
654
655Please report any bugs via:
656http://svn.atnf.csiro.au/trac/asap/simpleticket
657
658[IMPORTANT: ASAP is 0-based]
659Type commands() to get a list of all available ASAP commands.""" % (__version__, __date__)
Note: See TracBrowser for help on using the repository browser.