source: branches/Release-2-fixes/python/scantable.py@ 1898

Last change on this file since 1898 was 669, checked in by mar637, 19 years ago

fixed typo in import rotate_xyphase
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 57.2 KB
Line 
1from asap._asap import sdtable
2from asap import rcParams
3from numarray import ones,zeros
4import sys
5
6class scantable(sdtable):
7 """
8 The ASAP container for scans
9 """
10
11 def __init__(self, filename, average=None, unit=None):
12 """
13 Create a scantable from a saved one or make a reference
14 Parameters:
15 filename: the name of an asap table on disk
16 or
17 the name of a rpfits/sdfits/ms file
18 (integrations within scans are auto averaged
19 and the whole file is read)
20 or
21 [advanced] a reference to an existing
22 scantable
23 average: average all integrations withinb a scan on read.
24 The default (True) is taken from .asaprc.
25 unit: brightness unit; must be consistent with K or Jy.
26 Over-rides the default selected by the reader
27 (input rpfits/sdfits/ms) or replaces the value
28 in existing scantables
29 """
30 if average is None or type(average) is not bool:
31 average = rcParams['scantable.autoaverage']
32
33 varlist = vars()
34 self._vb = rcParams['verbose']
35 self._p = None
36
37 if isinstance(filename,sdtable):
38 sdtable.__init__(self, filename)
39 if unit is not None:
40 self.set_fluxunit(unit)
41 else:
42 import os.path
43 if not os.path.exists(filename):
44 print "File '%s' not found." % (filename)
45 return
46 filename = os.path.expandvars(filename)
47 if os.path.isdir(filename):
48 # crude check if asap table
49 if os.path.exists(filename+'/table.info'):
50 sdtable.__init__(self, filename)
51 if unit is not None:
52 self.set_fluxunit(unit)
53 else:
54 print "The given file '%s'is not a valid asap table." % (filename)
55 return
56 else:
57 from asap._asap import sdreader
58 ifSel = -1
59 beamSel = -1
60 r = sdreader(filename,ifSel,beamSel)
61 print 'Importing data...'
62 r._read([-1])
63 tbl = r._getdata()
64 if unit is not None:
65 tbl.set_fluxunit(unit)
66 if average:
67 from asap._asap import average as _av
68 print 'Auto averaging integrations...'
69 tbl2 = _av((tbl,),(),True,'none')
70 sdtable.__init__(self,tbl2)
71 del tbl2
72 else:
73 sdtable.__init__(self,tbl)
74 del r,tbl
75 self._add_history("scantable", varlist)
76
77 def save(self, name=None, format=None, stokes=False, overwrite=False):
78 """
79 Store the scantable on disk. This can be an asap (aips++) Table, SDFITS,
80 Image FITS or MS2 format.
81 Parameters:
82 name: the name of the outputfile. For format="FITS" this
83 is the directory file name into which all the files
84 will be written (default is 'asap_FITS'). For format
85 "ASCII" this is the root file name (data in 'name'.txt
86 and header in 'name'_header.txt)
87 format: an optional file format. Default is ASAP.
88 Allowed are - 'ASAP' (save as ASAP [aips++] Table),
89 'SDFITS' (save as SDFITS file)
90 'FITS' (saves each row as a FITS Image)
91 'ASCII' (saves as ascii text file)
92 'MS2' (saves as an aips++
93 MeasurementSet V2)
94 stokes: Convert to Stokes parameters (only available
95 currently with FITS and ASCII formats.
96 Default is False.
97 overwrite: If the file should be overwritten if it exists.
98 The default False is to return with warning
99 without writing the output. USE WITH CARE.
100 Example:
101 scan.save('myscan.asap')
102 scan.save('myscan.sdfits','SDFITS')
103 """
104 from os import path
105 if format is None: format = rcParams['scantable.save']
106 suffix = '.'+format.lower()
107 if name is None or name =="":
108 name = 'scantable'+suffix
109 print "No filename given. Using default name %s..." % name
110 name = path.expandvars(name)
111 if path.isfile(name) or path.isdir(name):
112 if not overwrite:
113 print "File %s already exists." % name
114 return
115 format2 = format.upper()
116 if format2 == 'ASAP':
117 self._save(name)
118 else:
119 from asap._asap import sdwriter as _sw
120 w = _sw(format2)
121 w.write(self, name, stokes)
122 return
123
124 def copy(self):
125 """
126 Return a copy of this scantable.
127 Parameters:
128 none
129 Example:
130 copiedscan = scan.copy()
131 """
132 sd = scantable(sdtable._copy(self))
133 return sd
134
135 def get_scan(self, scanid=None):
136 """
137 Return a specific scan (by scanno) or collection of scans (by
138 source name) in a new scantable.
139 Parameters:
140 scanid: a (list of) scanno or a source name, unix-style
141 patterns are accepted for source name matching, e.g.
142 '*_R' gets all 'ref scans
143 Example:
144 # get all scans containing the source '323p459'
145 newscan = scan.get_scan('323p459')
146 # get all 'off' scans
147 refscans = scan.get_scan('*_R')
148 # get a susbset of scans by scanno (as listed in scan.summary())
149 newscan = scan.get_scan([0,2,7,10])
150 """
151 if scanid is None:
152 print "Please specify a scan no or name to retrieve from the scantable"
153 try:
154 if type(scanid) is str:
155 s = sdtable._getsource(self,scanid)
156 return scantable(s)
157 elif type(scanid) is int:
158 s = sdtable._getscan(self,[scanid])
159 return scantable(s)
160 elif type(scanid) is list:
161 s = sdtable._getscan(self,scanid)
162 return scantable(s)
163 else:
164 print "Illegal scanid type, use 'int' or 'list' if ints."
165 except RuntimeError:
166 print "Couldn't find any match."
167
168 def __str__(self):
169 return sdtable._summary(self,True)
170
171 def summary(self,filename=None, verbose=None):
172 """
173 Print a summary of the contents of this scantable.
174 Parameters:
175 filename: the name of a file to write the putput to
176 Default - no file output
177 verbose: print extra info such as the frequency table
178 The default (False) is taken from .asaprc
179 """
180 info = sdtable._summary(self, verbose)
181 if verbose is None: verbose = rcParams['scantable.verbosesummary']
182 if filename is not None:
183 if filename is "":
184 filename = 'scantable_summary.txt'
185 from os.path import expandvars, isdir
186 filename = expandvars(filename)
187 if not isdir(filename):
188 data = open(filename, 'w')
189 data.write(info)
190 data.close()
191 else:
192 print "Illegal file name '%s'." % (filename)
193 print info
194
195 def set_cursor(self, beam=0, IF=0, pol=0):
196 """
197 Set the spectrum for individual operations.
198 Parameters:
199 beam, IF, pol: a number
200 Example:
201 scan.set_cursor(0,0,1)
202 pol1sig = scan.stats(all=False) # returns std dev for beam=0
203 # if=0, pol=1
204 """
205 varlist = vars()
206 self.setbeam(beam)
207 self.setpol(pol)
208 self.setif(IF)
209 self._add_history("set_cursor",varlist)
210 return
211
212 def get_cursor(self):
213 """
214 Return/print a the current 'cursor' into the Beam/IF/Pol cube.
215 Parameters:
216 none
217 Returns:
218 a list of values (currentBeam,currentIF,currentPol)
219 Example:
220 none
221 """
222 i = self.getbeam()
223 j = self.getif()
224 k = self.getpol()
225 if self._vb:
226 print "--------------------------------------------------"
227 print " Cursor position"
228 print "--------------------------------------------------"
229 out = 'Beam=%d IF=%d Pol=%d ' % (i,j,k)
230 print out
231 return i,j,k
232
233 def stats(self, stat='stddev', mask=None, allaxes=None):
234 """
235 Determine the specified statistic of the current beam/if/pol
236 Takes a 'mask' as an optional parameter to specify which
237 channels should be excluded.
238 Parameters:
239 stat: 'min', 'max', 'sumsq', 'sum', 'mean'
240 'var', 'stddev', 'avdev', 'rms', 'median'
241 mask: an optional mask specifying where the statistic
242 should be determined.
243 allaxes: if True apply to all spectra. Otherwise
244 apply only to the selected (beam/pol/if)spectra only.
245 The default is taken from .asaprc (True if none)
246 Example:
247 scan.set_unit('channel')
248 msk = scan.create_mask([100,200],[500,600])
249 scan.stats(stat='mean', mask=m)
250 """
251 if allaxes is None: allaxes = rcParams['scantable.allaxes']
252 from asap._asap import stats as _stats
253 from numarray import array,zeros,Float
254 if mask == None:
255 mask = ones(self.nchan())
256 axes = ['Beam','IF','Pol','Time']
257
258 beamSel,IFSel,polSel = (self.getbeam(),self.getif(),self.getpol())
259 if allaxes:
260 n = self.nbeam()*self.nif()*self.npol()*self.nrow()
261 shp = [self.nbeam(),self.nif(),self.npol(),self.nrow()]
262 arr = array(zeros(n),shape=shp,type=Float)
263
264 for i in range(self.nbeam()):
265 self.setbeam(i)
266 for j in range(self.nif()):
267 self.setif(j)
268 for k in range(self.npol()):
269 self.setpol(k)
270 arr[i,j,k,:] = _stats(self,mask,stat,-1)
271 retval = {'axes': axes, 'data': arr, 'cursor':None}
272 tm = [self._gettime(val) for val in range(self.nrow())]
273 if self._vb:
274 self._print_values(retval,stat,tm)
275 self.setbeam(beamSel)
276 self.setif(IFSel)
277 self.setpol(polSel)
278 return retval
279
280 else:
281 statval = _stats(self,mask,stat,-1)
282 out = ''
283 for l in range(self.nrow()):
284 tm = self._gettime(l)
285 out += 'Time[%s]:\n' % (tm)
286 if self.nbeam() > 1: out += ' Beam[%d] ' % (beamSel)
287 if self.nif() > 1: out += ' IF[%d] ' % (IFSel)
288 if self.npol() > 1: out += ' Pol[%d] ' % (polSel)
289 out += '= %3.3f\n' % (statval[l])
290 out += "--------------------------------------------------\n"
291
292 if self._vb:
293 print "--------------------------------------------------"
294 print " ",stat
295 print "--------------------------------------------------"
296 print out
297 retval = {'axes': axes, 'data': array(statval), 'cursor':(i,j,k)}
298 return retval
299
300 def stddev(self,mask=None, allaxes=None):
301 """
302 Determine the standard deviation of the current beam/if/pol
303 Takes a 'mask' as an optional parameter to specify which
304 channels should be excluded.
305 Parameters:
306 mask: an optional mask specifying where the standard
307 deviation should be determined.
308 allaxes: optional flag to show all or a cursor selected
309 spectrum of Beam/IF/Pol. Default is all or taken
310 from .asaprc
311
312 Example:
313 scan.set_unit('channel')
314 msk = scan.create_mask([100,200],[500,600])
315 scan.stddev(mask=m)
316 """
317 if allaxes is None: allaxes = rcParams['scantable.allaxes']
318 return self.stats(stat='stddev',mask=mask, allaxes=allaxes);
319
320 def get_tsys(self, allaxes=None):
321 """
322 Return the System temperatures.
323 Parameters:
324 allaxes: if True apply to all spectra. Otherwise
325 apply only to the selected (beam/pol/if)spectra only.
326 The default is taken from .asaprc (True if none)
327 Returns:
328 a list of Tsys values.
329 """
330 if allaxes is None: allaxes = rcParams['scantable.allaxes']
331 from numarray import array,zeros,Float
332 axes = ['Beam','IF','Pol','Time']
333
334 if allaxes:
335 n = self.nbeam()*self.nif()*self.npol()*self.nrow()
336 shp = [self.nbeam(),self.nif(),self.npol(),self.nrow()]
337 arr = array(zeros(n),shape=shp,type=Float)
338
339 for i in range(self.nbeam()):
340 self.setbeam(i)
341 for j in range(self.nif()):
342 self.setif(j)
343 for k in range(self.npol()):
344 self.setpol(k)
345 arr[i,j,k,:] = self._gettsys()
346 retval = {'axes': axes, 'data': arr, 'cursor':None}
347 tm = [self._gettime(val) for val in range(self.nrow())]
348 if self._vb:
349 self._print_values(retval,'Tsys',tm)
350 return retval
351
352 else:
353 i,j,k = (self.getbeam(),self.getif(),self.getpol())
354 statval = self._gettsys()
355 out = ''
356 for l in range(self.nrow()):
357 tm = self._gettime(l)
358 out += 'Time[%s]:\n' % (tm)
359 if self.nbeam() > 1: out += ' Beam[%d] ' % (i)
360 if self.nif() > 1: out += ' IF[%d] ' % (j)
361 if self.npol() > 1: out += ' Pol[%d] ' % (k)
362 out += '= %3.3f\n' % (statval[l])
363 out += "--------------------------------------------------\n"
364
365 if self._vb:
366 print "--------------------------------------------------"
367 print " TSys"
368 print "--------------------------------------------------"
369 print out
370 retval = {'axes': axes, 'data': array(statval), 'cursor':(i,j,k)}
371 return retval
372
373 def get_time(self, row=-1):
374 """
375 Get a list of time stamps for the observations.
376 Return a string for each integration in the scantable.
377 Parameters:
378 row: row no of integration. Default -1 return all rows
379 Example:
380 none
381 """
382 out = []
383 if row == -1:
384 for i in range(self.nrow()):
385 out.append(self._gettime(i))
386 return out
387 else:
388 if row < self.nrow():
389 return self._gettime(row)
390
391 def set_unit(self, unit='channel'):
392 """
393 Set the unit for all following operations on this scantable
394 Parameters:
395 unit: optional unit, default is 'channel'
396 one of '*Hz','km/s','channel', ''
397 """
398 varlist = vars()
399 if unit in ['','pixel', 'channel']:
400 unit = ''
401 inf = list(self._getcoordinfo())
402 inf[0] = unit
403 self._setcoordinfo(inf)
404 if self._p: self.plot()
405 self._add_history("set_unit",varlist)
406
407 def set_instrument(self, instr):
408 """
409 Set the instrument for subsequent processing
410 Parameters:
411 instr: Select from 'ATPKSMB', 'ATPKSHOH', 'ATMOPRA',
412 'DSS-43' (Tid), 'CEDUNA', and 'HOBART'
413 """
414 self._setInstrument(instr)
415 self._add_history("set_instument",vars())
416
417 def set_doppler(self, doppler='RADIO'):
418 """
419 Set the doppler for all following operations on this scantable.
420 Parameters:
421 doppler: One of 'RADIO', 'OPTICAL', 'Z', 'BETA', 'GAMMA'
422 """
423 varlist = vars()
424 inf = list(self._getcoordinfo())
425 inf[2] = doppler
426 self._setcoordinfo(inf)
427 if self._p: self.plot()
428 self._add_history("set_doppler",vars())
429
430 def set_freqframe(self, frame=None):
431 """
432 Set the frame type of the Spectral Axis.
433 Parameters:
434 frame: an optional frame type, default 'LSRK'. Valid frames are:
435 'REST','TOPO','LSRD','LSRK','BARY',
436 'GEO','GALACTO','LGROUP','CMB'
437 Examples:
438 scan.set_freqframe('BARY')
439 """
440 if frame is None: frame = rcParams['scantable.freqframe']
441 varlist = vars()
442 valid = ['REST','TOPO','LSRD','LSRK','BARY', \
443 'GEO','GALACTO','LGROUP','CMB']
444
445 if frame in valid:
446 inf = list(self._getcoordinfo())
447 inf[1] = frame
448 self._setcoordinfo(inf)
449 self._add_history("set_freqframe",varlist)
450 else:
451 print "Please specify a valid freq type. Valid types are:\n",valid
452
453 def get_unit(self):
454 """
455 Get the default unit set in this scantable
456 Parameters:
457 Returns:
458 A unit string
459 """
460 inf = self._getcoordinfo()
461 unit = inf[0]
462 if unit == '': unit = 'channel'
463 return unit
464
465 def get_abcissa(self, rowno=0):
466 """
467 Get the abcissa in the current coordinate setup for the currently
468 selected Beam/IF/Pol
469 Parameters:
470 rowno: an optional row number in the scantable. Default is the
471 first row, i.e. rowno=0
472 Returns:
473 The abcissa values and it's format string (as a dictionary)
474 """
475 abc = self._getabcissa(rowno)
476 lbl = self._getabcissalabel(rowno)
477 return abc, lbl
478 #return {'abcissa':abc,'label':lbl}
479
480 def create_mask(self, *args, **kwargs):
481 """
482 Compute and return a mask based on [min,max] windows.
483 The specified windows are to be INCLUDED, when the mask is
484 applied.
485 Parameters:
486 [min,max],[min2,max2],...
487 Pairs of start/end points specifying the regions
488 to be masked
489 invert: optional argument. If specified as True,
490 return an inverted mask, i.e. the regions
491 specified are EXCLUDED
492 row: create the mask using the specified row for
493 unit conversions, default is row=0
494 only necessary if frequency varies over rows.
495 Example:
496 scan.set_unit('channel')
497
498 a)
499 msk = scan.set_mask([400,500],[800,900])
500 # masks everything outside 400 and 500
501 # and 800 and 900 in the unit 'channel'
502
503 b)
504 msk = scan.set_mask([400,500],[800,900], invert=True)
505 # masks the regions between 400 and 500
506 # and 800 and 900 in the unit 'channel'
507
508 """
509 row = 0
510 if kwargs.has_key("row"):
511 row = kwargs.get("row")
512 data = self._getabcissa(row)
513 u = self._getcoordinfo()[0]
514 if self._vb:
515 if u == "": u = "channel"
516 print "The current mask window unit is", u
517 n = self.nchan()
518 msk = zeros(n)
519 for window in args:
520 if (len(window) != 2 or window[0] > window[1] ):
521 print "A window needs to be defined as [min,max]"
522 return
523 for i in range(n):
524 if data[i] >= window[0] and data[i] < window[1]:
525 msk[i] = 1
526 if kwargs.has_key('invert'):
527 if kwargs.get('invert'):
528 from numarray import logical_not
529 msk = logical_not(msk)
530 return msk
531
532 def get_restfreqs(self):
533 """
534 Get the restfrequency(s) stored in this scantable.
535 The return value(s) are always of unit 'Hz'
536 Parameters:
537 none
538 Returns:
539 a list of doubles
540 """
541 return list(self._getrestfreqs())
542
543 def lines(self):
544 """
545 Print the list of known spectral lines
546 """
547 sdtable._lines(self)
548
549 def set_restfreqs(self, freqs=None, unit='Hz', lines=None, source=None,
550 theif=None):
551 """
552 Select the restfrequency for the specified source and IF OR
553 replace for all IFs. If the 'freqs' argument holds a scalar,
554 then that rest frequency will be applied to the selected
555 data (and added to the list of available rest frequencies).
556 In this way, you can set a rest frequency for each
557 source and IF combination. If the 'freqs' argument holds
558 a vector, then it MUST be of length the number of IFs
559 (and the available restfrequencies will be replaced by
560 this vector). In this case, *all* data ('source' and
561 'theif' are ignored) have the restfrequency set per IF according
562 to the corresponding value you give in the 'freqs' vector.
563 E.g. 'freqs=[1e9,2e9]' would mean IF 0 gets restfreq 1e9 and
564 IF 1 gets restfreq 2e9.
565
566 You can also specify the frequencies via known line names
567 in the argument 'lines'. Use 'freqs' or 'lines'. 'freqs'
568 takes precedence. See the list of known names via function
569 scantable.lines()
570 Parameters:
571 freqs: list of rest frequencies
572 unit: unit for rest frequency (default 'Hz')
573 lines: list of known spectral lines names (alternative to freqs).
574 See possible list via scantable.lines()
575 source: Source name (blank means all)
576 theif: IF (-1 means all)
577 Example:
578 scan.set_restfreqs(freqs=1.4e9, source='NGC253', theif=2)
579 scan.set_restfreqs(freqs=[1.4e9,1.67e9])
580 """
581 varlist = vars()
582 if source is None:
583 source = ""
584 if theif is None:
585 theif = -1
586 t = type(freqs)
587 if t is int or t is float:
588 freqs = [freqs]
589 if freqs is None:
590 freqs = []
591 t = type(lines)
592 if t is str:
593 lines = [lines]
594 if lines is None:
595 lines = []
596 sdtable._setrestfreqs(self, freqs, unit, lines, source, theif)
597 self._add_history("set_restfreqs", varlist)
598
599
600
601 def flag_spectrum(self, thebeam, theif, thepol):
602 """
603 This flags a selected spectrum in the scan 'for good'.
604 USE WITH CARE - not reversible.
605 Use masks for non-permanent exclusion of channels.
606 Parameters:
607 thebeam,theif,thepol: all have to be explicitly
608 specified
609 Example:
610 scan.flag_spectrum(0,0,1)
611 flags the spectrum for Beam=0, IF=0, Pol=1
612 """
613 if (thebeam < self.nbeam() and
614 theif < self.nif() and
615 thepol < self.npol()):
616 sdtable.setbeam(self, thebeam)
617 sdtable.setif(self, theif)
618 sdtable.setpol(self, thepol)
619 sdtable._flag(self)
620 self._add_history("flag_spectrum", vars())
621 else:
622 print "Please specify a valid (Beam/IF/Pol)"
623 return
624
625 def plot(self, what='spectrum',col='Pol', panel=None):
626 """
627 Plot the spectra contained in the scan. Alternatively you can also
628 Plot Tsys vs Time
629 Parameters:
630 what: a choice of 'spectrum' (default) or 'tsys'
631 col: which out of Beams/IFs/Pols should be colour stacked
632 panel: set up multiple panels, currently not working.
633 """
634 print "Warning! Not fully functional. Use plotter.plot() instead"
635
636 validcol = {'Beam':self.nbeam(),'IF':self.nif(),'Pol':self.npol()}
637
638 validyax = ['spectrum','tsys']
639 from asap.asaplot import ASAPlot
640 if not self._p:
641 self._p = ASAPlot()
642 #print "Plotting not enabled"
643 #return
644 if self._p.is_dead:
645 del self._p
646 self._p = ASAPlot()
647 npan = 1
648 x = None
649 if what == 'tsys':
650 n = self.nrow()
651 if n < 2:
652 print "Only one integration. Can't plot."
653 return
654 self._p.hold()
655 self._p.clear()
656 if panel == 'Time':
657 npan = self.nrow()
658 self._p.set_panels(rows=npan)
659 xlab,ylab,tlab = None,None,None
660 self._vb = False
661 sel = self.get_cursor()
662 for i in range(npan):
663 if npan > 1:
664 self._p.subplot(i)
665 for j in range(validcol[col]):
666 x = None
667 y = None
668 m = None
669 tlab = self._getsourcename(i)
670 import re
671 tlab = re.sub('_S','',tlab)
672 if col == 'Beam':
673 self.setbeam(j)
674 elif col == 'IF':
675 self.setif(j)
676 elif col == 'Pol':
677 self.setpol(j)
678 if what == 'tsys':
679 x = range(self.nrow())
680 xlab = 'Time [pixel]'
681 m = list(ones(len(x)))
682 y = []
683 ylab = r'$T_{sys}$'
684 for k in range(len(x)):
685 y.append(self._gettsys(k))
686 else:
687 x,xlab = self.get_abcissa(i)
688 y = self._getspectrum(i)
689 ylab = r'Flux'
690 m = self._getmask(i)
691 llab = col+' '+str(j)
692 self._p.set_line(label=llab)
693 self._p.plot(x,y,m)
694 self._p.set_axes('xlabel',xlab)
695 self._p.set_axes('ylabel',ylab)
696 self._p.set_axes('title',tlab)
697 self._p.release()
698 self.set_cursor(sel[0],sel[1],sel[2])
699 self._vb = rcParams['verbose']
700 return
701
702 print out
703
704 def _print_values(self, dat, label='', timestamps=[]):
705 d = dat['data']
706 a = dat['axes']
707 shp = d.getshape()
708 out = ''
709 for i in range(shp[3]):
710 out += '%s [%s]:\n' % (a[3],timestamps[i])
711 t = d[:,:,:,i]
712 for j in range(shp[0]):
713 if shp[0] > 1: out += ' %s[%d] ' % (a[0],j)
714 for k in range(shp[1]):
715 if shp[1] > 1: out += ' %s[%d] ' % (a[1],k)
716 for l in range(shp[2]):
717 if shp[2] > 1: out += ' %s[%d] ' % (a[2],l)
718 out += '= %3.3f\n' % (t[j,k,l])
719 out += "-"*80
720 out += "\n"
721 print "-"*80
722 print " ", label
723 print "-"*80
724 print out
725
726 def history(self):
727 hist = list(self._gethistory())
728 print "-"*80
729 for h in hist:
730 if h.startswith("---"):
731 print h
732 else:
733 items = h.split("##")
734 date = items[0]
735 func = items[1]
736 items = items[2:]
737 print date
738 print "Function: %s\n Parameters:" % (func)
739 for i in items:
740 s = i.split("=")
741 print " %s = %s" % (s[0],s[1])
742 print "-"*80
743 return
744
745 #
746 # Maths business
747 #
748
749 def average_time(self, mask=None, scanav=False, weight='tint'):
750 """
751 Return the (time) average of a scan, or apply it 'insitu'.
752 Note:
753 in channels only
754 The cursor of the output scan is set to 0.
755 Parameters:
756 one scan or comma separated scans
757 mask: an optional mask (only used for 'var' and 'tsys'
758 weighting)
759 scanav: True averages each scan separately
760 False (default) averages all scans together,
761 weight: Weighting scheme. 'none', 'var' (1/var(spec)
762 weighted), 'tsys' (1/Tsys**2 weighted), 'tint'
763 (integration time weighted) or 'tintsys' (Tint/Tsys**2).
764 The default is 'tint'
765 Example:
766 # time average the scantable without using a mask
767 newscan = scan.average_time()
768 """
769 varlist = vars()
770 if weight is None: weight = 'tint'
771 if mask is None: mask = ()
772 from asap._asap import average as _av
773 s = scantable(_av((self,), mask, scanav, weight))
774 s._add_history("average_time",varlist)
775 return s
776
777 def convert_flux(self, jyperk=None, eta=None, d=None, insitu=None,
778 allaxes=None):
779 """
780 Return a scan where all spectra are converted to either
781 Jansky or Kelvin depending upon the flux units of the scan table.
782 By default the function tries to look the values up internally.
783 If it can't find them (or if you want to over-ride), you must
784 specify EITHER jyperk OR eta (and D which it will try to look up
785 also if you don't set it). jyperk takes precedence if you set both.
786 Parameters:
787 jyperk: the Jy / K conversion factor
788 eta: the aperture efficiency
789 d: the geomtric diameter (metres)
790 insitu: if False a new scantable is returned.
791 Otherwise, the scaling is done in-situ
792 The default is taken from .asaprc (False)
793 allaxes: if True apply to all spectra. Otherwise
794 apply only to the selected (beam/pol/if)spectra only
795 The default is taken from .asaprc (True if none)
796 """
797 if allaxes is None: allaxes = rcParams['scantable.allaxes']
798 if insitu is None: insitu = rcParams['insitu']
799 varlist = vars()
800 if jyperk is None: jyperk = -1.0
801 if d is None: d = -1.0
802 if eta is None: eta = -1.0
803 if not insitu:
804 from asap._asap import convertflux as _convert
805 s = scantable(_convert(self, d, eta, jyperk, allaxes))
806 s._add_history("convert_flux", varlist)
807 return s
808 else:
809 from asap._asap import convertflux_insitu as _convert
810 _convert(self, d, eta, jyperk, allaxes)
811 self._add_history("convert_flux", varlist)
812 return
813
814 def gain_el(self, poly=None, filename="", method="linear",
815 insitu=None, allaxes=None):
816 """
817 Return a scan after applying a gain-elevation correction.
818 The correction can be made via either a polynomial or a
819 table-based interpolation (and extrapolation if necessary).
820 You specify polynomial coefficients, an ascii table or neither.
821 If you specify neither, then a polynomial correction will be made
822 with built in coefficients known for certain telescopes (an error
823 will occur if the instrument is not known).
824 The data and Tsys are *divided* by the scaling factors.
825 Parameters:
826 poly: Polynomial coefficients (default None) to compute a
827 gain-elevation correction as a function of
828 elevation (in degrees).
829 filename: The name of an ascii file holding correction factors.
830 The first row of the ascii file must give the column
831 names and these MUST include columns
832 "ELEVATION" (degrees) and "FACTOR" (multiply data
833 by this) somewhere.
834 The second row must give the data type of the
835 column. Use 'R' for Real and 'I' for Integer.
836 An example file would be
837 (actual factors are arbitrary) :
838
839 TIME ELEVATION FACTOR
840 R R R
841 0.1 0 0.8
842 0.2 20 0.85
843 0.3 40 0.9
844 0.4 60 0.85
845 0.5 80 0.8
846 0.6 90 0.75
847 method: Interpolation method when correcting from a table.
848 Values are "nearest", "linear" (default), "cubic"
849 and "spline"
850 insitu: if False a new scantable is returned.
851 Otherwise, the scaling is done in-situ
852 The default is taken from .asaprc (False)
853 allaxes: If True apply to all spectra. Otherwise
854 apply only to the selected (beam/pol/if) spectra only
855 The default is taken from .asaprc (True if none)
856 """
857
858 if allaxes is None: allaxes = rcParams['scantable.allaxes']
859 if insitu is None: insitu = rcParams['insitu']
860 varlist = vars()
861 if poly is None:
862 poly = ()
863 from os.path import expandvars
864 filename = expandvars(filename)
865 if not insitu:
866 from asap._asap import gainel as _gainEl
867 s = scantable(_gainEl(self, poly, filename, method, allaxes))
868 s._add_history("gain_el", varlist)
869 return s
870 else:
871 from asap._asap import gainel_insitu as _gainEl
872 _gainEl(self, poly, filename, method, allaxes)
873 self._add_history("gain_el", varlist)
874 return
875
876 def freq_align(self, reftime=None, method='cubic', perif=False,
877 insitu=None):
878 """
879 Return a scan where all rows have been aligned in frequency/velocity.
880 The alignment frequency frame (e.g. LSRK) is that set by function
881 set_freqframe.
882 Parameters:
883 reftime: reference time to align at. By default, the time of
884 the first row of data is used.
885 method: Interpolation method for regridding the spectra.
886 Choose from "nearest", "linear", "cubic" (default)
887 and "spline"
888 perif: Generate aligners per freqID (no doppler tracking) or
889 per IF (scan-based doppler tracking)
890 insitu: if False a new scantable is returned.
891 Otherwise, the scaling is done in-situ
892 The default is taken from .asaprc (False)
893 """
894 if insitu is None: insitu = rcParams['insitu']
895 varlist = vars()
896 if reftime is None: reftime = ''
897 perfreqid = not perif
898 if not insitu:
899 from asap._asap import freq_align as _align
900 s = scantable(_align(self, reftime, method, perfreqid))
901 s._add_history("freq_align", varlist)
902 return s
903 else:
904 from asap._asap import freq_align_insitu as _align
905 _align(self, reftime, method, perfreqid)
906 self._add_history("freq_align", varlist)
907 return
908
909 def opacity(self, tau, insitu=None, allaxes=None):
910 """
911 Apply an opacity correction. The data
912 and Tsys are multiplied by the correction factor.
913 Parameters:
914 tau: Opacity from which the correction factor is
915 exp(tau*ZD)
916 where ZD is the zenith-distance
917 insitu: if False a new scantable is returned.
918 Otherwise, the scaling is done in-situ
919 The default is taken from .asaprc (False)
920 allaxes: if True apply to all spectra. Otherwise
921 apply only to the selected (beam/pol/if)spectra only
922 The default is taken from .asaprc (True if none)
923 """
924 if allaxes is None: allaxes = rcParams['scantable.allaxes']
925 if insitu is None: insitu = rcParams['insitu']
926 varlist = vars()
927 if not insitu:
928 from asap._asap import opacity as _opacity
929 s = scantable(_opacity(self, tau, allaxes))
930 s._add_history("opacity", varlist)
931 return s
932 else:
933 from asap._asap import opacity_insitu as _opacity
934 _opacity(self, tau, allaxes)
935 self._add_history("opacity", varlist)
936 return
937
938 def bin(self, width=5, insitu=None):
939 """
940 Return a scan where all spectra have been binned up.
941 width: The bin width (default=5) in pixels
942 insitu: if False a new scantable is returned.
943 Otherwise, the scaling is done in-situ
944 The default is taken from .asaprc (False)
945 """
946 if insitu is None: insitu = rcParams['insitu']
947 varlist = vars()
948 if not insitu:
949 from asap._asap import bin as _bin
950 s = scantable(_bin(self, width))
951 s._add_history("bin",varlist)
952 return s
953 else:
954 from asap._asap import bin_insitu as _bin
955 _bin(self, width)
956 self._add_history("bin",varlist)
957 return
958
959
960 def resample(self, width=5, method='cubic', insitu=None):
961 """
962 Return a scan where all spectra have been binned up
963 width: The bin width (default=5) in pixels
964 method: Interpolation method when correcting from a table.
965 Values are "nearest", "linear", "cubic" (default)
966 and "spline"
967 insitu: if False a new scantable is returned.
968 Otherwise, the scaling is done in-situ
969 The default is taken from .asaprc (False)
970 """
971 if insitu is None: insitu = rcParams['insitu']
972 varlist = vars()
973 if not insitu:
974 from asap._asap import resample as _resample
975 s = scantable(_resample(self, method, width))
976 s._add_history("resample",varlist)
977 return s
978 else:
979 from asap._asap import resample_insitu as _resample
980 _resample(self, method, width)
981 self._add_history("resample",varlist)
982 return
983
984 def average_pol(self, mask=None, weight='none', insitu=None):
985 """
986 Average the Polarisations together.
987 The polarisation cursor of the output scan is set to 0
988 Parameters:
989 mask: An optional mask defining the region, where the
990 averaging will be applied. The output will have all
991 specified points masked.
992 weight: Weighting scheme. 'none' (default), 'var' (1/var(spec)
993 weighted), or 'tsys' (1/Tsys**2 weighted)
994 insitu: if False a new scantable is returned.
995 Otherwise, the scaling is done in-situ
996 The default is taken from .asaprc (False)
997 """
998 if insitu is None: insitu = rcParams['insitu']
999 varlist = vars()
1000
1001 if mask is None:
1002 mask = ()
1003 if not insitu:
1004 from asap._asap import averagepol as _avpol
1005 s = scantable(_avpol(self, mask, weight))
1006 s._add_history("average_pol",varlist)
1007 return s
1008 else:
1009 from asap._asap import averagepol_insitu as _avpol
1010 _avpol(self, mask, weight)
1011 self._add_history("average_pol",varlist)
1012 return
1013
1014 def smooth(self, kernel="hanning", width=5.0, insitu=None, allaxes=None):
1015 """
1016 Smooth the spectrum by the specified kernel (conserving flux).
1017 Parameters:
1018 scan: The input scan
1019 kernel: The type of smoothing kernel. Select from
1020 'hanning' (default), 'gaussian' and 'boxcar'.
1021 The first three characters are sufficient.
1022 width: The width of the kernel in pixels. For hanning this is
1023 ignored otherwise it defauls to 5 pixels.
1024 For 'gaussian' it is the Full Width Half
1025 Maximum. For 'boxcar' it is the full width.
1026 insitu: if False a new scantable is returned.
1027 Otherwise, the scaling is done in-situ
1028 The default is taken from .asaprc (False)
1029 allaxes: If True (default) apply to all spectra. Otherwise
1030 apply only to the selected (beam/pol/if)spectra only
1031 The default is taken from .asaprc (True if none)
1032 Example:
1033 none
1034 """
1035 if allaxes is None: allaxes = rcParams['scantable.allaxes']
1036 if insitu is None: insitu = rcParams['insitu']
1037 varlist = vars()
1038 if not insitu:
1039 from asap._asap import smooth as _smooth
1040 s = scantable(_smooth(self,kernel,width,allaxes))
1041 s._add_history("smooth", varlist)
1042 return s
1043 else:
1044 from asap._asap import smooth_insitu as _smooth
1045 _smooth(self,kernel,width,allaxes)
1046 self._add_history("smooth", varlist)
1047 return
1048
1049 def poly_baseline(self, mask=None, order=0, insitu=None):
1050 """
1051 Return a scan which has been baselined (all rows) by a polynomial.
1052 Parameters:
1053 scan: a scantable
1054 mask: an optional mask
1055 order: the order of the polynomial (default is 0)
1056 insitu: if False a new scantable is returned.
1057 Otherwise, the scaling is done in-situ
1058 The default is taken from .asaprc (False)
1059 Example:
1060 # return a scan baselined by a third order polynomial,
1061 # not using a mask
1062 bscan = scan.poly_baseline(order=3)
1063 """
1064 if insitu is None: insitu = rcParams['insitu']
1065 varlist = vars()
1066 if mask is None:
1067 from numarray import ones
1068 mask = list(ones(self.nchan()))
1069 from asap.asapfitter import fitter
1070 f = fitter()
1071 f._verbose(True)
1072 f.set_scan(self, mask)
1073 f.set_function(poly=order)
1074 sf = f.auto_fit(insitu)
1075 if insitu:
1076 self._add_history("poly_baseline", varlist)
1077 return
1078 else:
1079 sf._add_history("poly_baseline", varlist)
1080 return sf
1081
1082 def auto_poly_baseline(self, mask=None, edge=(0,0), order=0,
1083 threshold=3,insitu=None):
1084 """
1085 Return a scan which has been baselined (all rows) by a polynomial.
1086 Spectral lines are detected first using linefinder and masked out
1087 to avoid them affecting the baseline solution.
1088
1089 Parameters:
1090 scan: a scantable
1091 mask: an optional mask retreived from scantable
1092 edge: an optional number of channel to drop at
1093 the edge of spectrum. If only one value is
1094 specified, the same number will be dropped from
1095 both sides of the spectrum. Default is to keep
1096 all channels
1097 order: the order of the polynomial (default is 0)
1098 threshold: the threshold used by line finder. It is better to
1099 keep it large as only strong lines affect the
1100 baseline solution.
1101 insitu: if False a new scantable is returned.
1102 Otherwise, the scaling is done in-situ
1103 The default is taken from .asaprc (False)
1104
1105 Example:
1106 scan2=scan.auto_poly_baseline(order=7)
1107 """
1108 if insitu is None: insitu = rcParams['insitu']
1109 varlist = vars()
1110 from asap.asapfitter import fitter
1111 from asap.asaplinefind import linefinder
1112 from asap import _is_sequence_or_number as _is_valid
1113
1114 if not _is_valid(edge, int):
1115 raise RuntimeError, "Parameter 'edge' has to be an integer or a \
1116 pair of integers specified as a tuple"
1117
1118 # setup fitter
1119 f = fitter()
1120 f._verbose(True)
1121 f.set_function(poly=order)
1122
1123 # setup line finder
1124 fl=linefinder()
1125 fl.set_options(threshold=threshold)
1126
1127 if not insitu:
1128 workscan=self.copy()
1129 else:
1130 workscan=self
1131
1132 vb=workscan._vb
1133 # remember the verbose parameter and selection
1134 workscan._vb=False
1135 sel=workscan.get_cursor()
1136 rows=range(workscan.nrow())
1137 for i in range(workscan.nbeam()):
1138 workscan.setbeam(i)
1139 for j in range(workscan.nif()):
1140 workscan.setif(j)
1141 for k in range(workscan.npol()):
1142 workscan.setpol(k)
1143 if f._vb:
1144 print "Processing:"
1145 print 'Beam[%d], IF[%d], Pol[%d]' % (i,j,k)
1146 for iRow in rows:
1147 fl.set_scan(workscan,mask,edge)
1148 fl.find_lines(iRow)
1149 f.set_scan(workscan, fl.get_mask())
1150 f.x=workscan._getabcissa(iRow)
1151 f.y=workscan._getspectrum(iRow)
1152 f.data=None
1153 f.fit()
1154 x=f.get_parameters()
1155 workscan._setspectrum(f.fitter.getresidual(),iRow)
1156 workscan.set_cursor(sel[0],sel[1],sel[2])
1157 workscan._vb = vb
1158 if not insitu:
1159 return workscan
1160
1161 def rotate_linpolphase(self, angle, allaxes=None):
1162 """
1163 Rotate the phase of the complex polarization O=Q+iU correlation.
1164 This is always done in situ in the raw data. So if you call this
1165 function more than once then each call rotates the phase further.
1166 Parameters:
1167 angle: The angle (degrees) to rotate (add) by.
1168 allaxes: If True apply to all spectra. Otherwise
1169 apply only to the selected (beam/pol/if)spectra only.
1170 The default is taken from .asaprc (True if none)
1171 Examples:
1172 scan.rotate_linpolphase(2.3)
1173 """
1174 if allaxes is None: allaxes = rcParams['scantable.allaxes']
1175 varlist = vars()
1176 from asap._asap import _rotate_linpolphase as _rotate
1177 _rotate(self, angle, allaxes)
1178 self._add_history("rotate_linpolphase", varlist)
1179 return
1180
1181
1182 def rotate_xyphase(self, angle, allaxes=None):
1183 """
1184 Rotate the phase of the XY correlation. This is always done in situ
1185 in the data. So if you call this function more than once
1186 then each call rotates the phase further.
1187 Parameters:
1188 angle: The angle (degrees) to rotate (add) by.
1189 allaxes: If True apply to all spectra. Otherwise
1190 apply only to the selected (beam/pol/if)spectra only.
1191 The default is taken from .asaprc (True if none)
1192 Examples:
1193 scan.rotate_xyphase(2.3)
1194 """
1195 if allaxes is None: allaxes = rcParams['scantable.allaxes']
1196 varlist = vars()
1197 from asap._asap import _rotate_xyphase
1198 _rotate_xyphase(self, angle, allaxes)
1199 self._add_history("rotate_xyphase", varlist)
1200 return
1201
1202
1203 def add(self, offset, insitu=None, allaxes=None):
1204 """
1205 Return a scan where all spectra have the offset added
1206 Parameters:
1207 offset: the offset
1208 insitu: if False a new scantable is returned.
1209 Otherwise, the scaling is done in-situ
1210 The default is taken from .asaprc (False)
1211 allaxes: if True apply to all spectra. Otherwise
1212 apply only to the selected (beam/pol/if)spectra only
1213 The default is taken from .asaprc (True if none)
1214 """
1215 if allaxes is None: allaxes = rcParams['scantable.allaxes']
1216 if insitu is None: insitu = rcParams['insitu']
1217 varlist = vars()
1218 if not insitu:
1219 from asap._asap import add as _add
1220 s = scantable(_add(self, offset, allaxes))
1221 s._add_history("add",varlist)
1222 return s
1223 else:
1224 from asap._asap import add_insitu as _add
1225 _add(self, offset, allaxes)
1226 self._add_history("add",varlist)
1227 return
1228
1229 def scale(self, factor, insitu=None, allaxes=None, tsys=True):
1230 """
1231 Return a scan where all spectra are scaled by the give 'factor'
1232 Parameters:
1233 factor: the scaling factor
1234 insitu: if False a new scantable is returned.
1235 Otherwise, the scaling is done in-situ
1236 The default is taken from .asaprc (False)
1237 allaxes: if True apply to all spectra. Otherwise
1238 apply only to the selected (beam/pol/if)spectra only.
1239 The default is taken from .asaprc (True if none)
1240 tsys: if True (default) then apply the operation to Tsys
1241 as well as the data
1242 """
1243 if allaxes is None: allaxes = rcParams['scantable.allaxes']
1244 if insitu is None: insitu = rcParams['insitu']
1245 varlist = vars()
1246 if not insitu:
1247 from asap._asap import scale as _scale
1248 s = scantable(_scale(self, factor, allaxes, tsys))
1249 s._add_history("scale",varlist)
1250 return s
1251 else:
1252 from asap._asap import scale_insitu as _scale
1253 _scale(self, factor, allaxes, tsys)
1254 self._add_history("scale",varlist)
1255 return
1256
1257 def auto_quotient(self, mode='suffix', preserve=True):
1258 """
1259 This function allows to build quotients automatically.
1260 It assumes the observation to have the same numer of
1261 "ons" and "offs"
1262 It will support "closest off in time" in the future
1263 Parameters:
1264 mode: the on/off detection mode; 'suffix' (default)
1265 'suffix' identifies 'off' scans by the
1266 trailing '_R' (Mopra/Parkes) or
1267 '_e'/'_w' (Tid)
1268 preserve: you can preserve (default) the continuum or
1269 remove it. The equations used are
1270 preserve: Output = Toff * (on/off) - Toff
1271 remove: Output = Tref * (on/off) - Ton
1272 """
1273 modes = ["suffix","time"]
1274 print mode
1275 if not mode in modes:
1276 print "please provide valid mode. Valid modes are %s" % (modes)
1277 return None
1278 from asap._asap import quotient as _quot
1279 if mode == "suffix":
1280 srcs = self.get_scan("*[^_ewR]")
1281 refs = self.get_scan("*[_ewR]")
1282 return scantable(_quot(srcs,refs, preserve))
1283 else:
1284 print "not yet implemented"
1285 return None
1286
1287 def quotient(self, other, isreference=True, preserve=True):
1288 """
1289 Return the quotient of a 'source' (on) scan and a 'reference' (off)
1290 scan.
1291 The reference can have just one row, even if the signal has many.
1292 Otherwise they must have the same number of rows.
1293 The cursor of the output scan is set to 0
1294 Parameters:
1295 other: the 'other' scan
1296 isreference: if the 'other' scan is the reference (default)
1297 or source
1298 preserve: you can preserve (default) the continuum or
1299 remove it. The equations used are
1300 preserve: Output = Toff * (on/off) - Toff
1301 remove: Output = Tref * (on/off) - Ton
1302 Example:
1303 # src is a scantable for an 'on' scan, ref for an 'off' scantable
1304 q1 = src.quotient(ref)
1305 q2 = ref.quotient(src, isreference=False)
1306 # gives the same result
1307 """
1308 from asap._asap import quotient as _quot
1309 if isreference:
1310 return scantable(_quot(self, other, preserve))
1311 else:
1312 return scantable(_quot(other, self, preserve))
1313
1314 def __add__(self, other):
1315 varlist = vars()
1316 s = None
1317 if isinstance(other, scantable):
1318 from asap._asap import b_operate as _bop
1319 s = scantable(_bop(self, other, 'add', True))
1320 elif isinstance(other, float):
1321 from asap._asap import add as _add
1322 s = scantable(_add(self, other, True))
1323 else:
1324 print "Other input is not a scantable or float value"
1325 return
1326 s._add_history("operator +", varlist)
1327 return s
1328
1329
1330 def __sub__(self, other):
1331 """
1332 implicit on all axes and on Tsys
1333 """
1334 varlist = vars()
1335 s = None
1336 if isinstance(other, scantable):
1337 from asap._asap import b_operate as _bop
1338 s = scantable(_bop(self, other, 'sub', True))
1339 elif isinstance(other, float):
1340 from asap._asap import add as _add
1341 s = scantable(_add(self, -other, True))
1342 else:
1343 print "Other input is not a scantable or float value"
1344 return
1345 s._add_history("operator -", varlist)
1346 return s
1347
1348 def __mul__(self, other):
1349 """
1350 implicit on all axes and on Tsys
1351 """
1352 varlist = vars()
1353 s = None
1354 if isinstance(other, scantable):
1355 from asap._asap import b_operate as _bop
1356 s = scantable(_bop(self, other, 'mul', True))
1357 elif isinstance(other, float):
1358 if other == 0.0:
1359 print "Multiplying by zero is not recommended."
1360 return
1361 from asap._asap import scale as _sca
1362 s = scantable(_sca(self, other, True))
1363 else:
1364 print "Other input is not a scantable or float value"
1365 return
1366 s._add_history("operator *", varlist)
1367 return s
1368
1369
1370 def __div__(self, other):
1371 """
1372 implicit on all axes and on Tsys
1373 """
1374 varlist = vars()
1375 s = None
1376 if isinstance(other, scantable):
1377 from asap._asap import b_operate as _bop
1378 s = scantable(_bop(self, other, 'div', True))
1379 elif isinstance(other, float):
1380 if other == 0.0:
1381 print "Dividing by zero is not recommended"
1382 return
1383 from asap._asap import scale as _sca
1384 s = scantable(_sca(self, 1.0/other, True))
1385 else:
1386 print "Other input is not a scantable or float value"
1387 return
1388 s._add_history("operator /", varlist)
1389 return s
1390
1391 def get_fit(self, row=0):
1392 """
1393 Print or return the stored fits for a row in the scantable
1394 Parameters:
1395 row: the row which the fit has been applied to.
1396 """
1397 if row > self.nrow():
1398 return
1399 from asap import asapfit
1400 fit = asapfit(self._getfit(row))
1401 if self._vb:
1402 print fit
1403 return
1404 else:
1405 return fit.as_dict()
1406
1407 def _add_history(self, funcname, parameters):
1408 # create date
1409 sep = "##"
1410 from datetime import datetime
1411 dstr = datetime.now().strftime('%Y/%m/%d %H:%M:%S')
1412 hist = dstr+sep
1413 hist += funcname+sep#cdate+sep
1414 if parameters.has_key('self'): del parameters['self']
1415 for k,v in parameters.iteritems():
1416 if type(v) is dict:
1417 for k2,v2 in v.iteritems():
1418 hist += k2
1419 hist += "="
1420 if isinstance(v2,scantable):
1421 hist += 'scantable'
1422 elif k2 == 'mask':
1423 if isinstance(v2,list) or isinstance(v2,tuple):
1424 hist += str(self._zip_mask(v2))
1425 else:
1426 hist += str(v2)
1427 else:
1428 hist += str(v2)
1429 else:
1430 hist += k
1431 hist += "="
1432 if isinstance(v,scantable):
1433 hist += 'scantable'
1434 elif k == 'mask':
1435 if isinstance(v,list) or isinstance(v,tuple):
1436 hist += str(self._zip_mask(v))
1437 else:
1438 hist += str(v)
1439 else:
1440 hist += str(v)
1441 hist += sep
1442 hist = hist[:-2] # remove trailing '##'
1443 self._addhistory(hist)
1444
1445
1446 def _zip_mask(self, mask):
1447 mask = list(mask)
1448 i = 0
1449 segments = []
1450 while mask[i:].count(1):
1451 i += mask[i:].index(1)
1452 if mask[i:].count(0):
1453 j = i + mask[i:].index(0)
1454 else:
1455 j = len(mask)
1456 segments.append([i,j])
1457 i = j
1458 return segments
1459 def _get_ordinate_label(self):
1460 fu = "("+self.get_fluxunit()+")"
1461 import re
1462 lbl = "Intensity"
1463 if re.match(".K.",fu):
1464 lbl = "Brightness Temperature "+ fu
1465 elif re.match(".Jy.",fu):
1466 lbl = "Flux density "+ fu
1467 return lbl
1468
Note: See TracBrowser for help on using the repository browser.