source: trunk/python/scantable.py@ 644

Last change on this file since 644 was 626, checked in by mar637, 20 years ago

fix for asap0019 from Release-1

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 55.9 KB
RevLine 
[102]1from asap._asap import sdtable
[226]2from asap import rcParams
[189]3from numarray import ones,zeros
[102]4import sys
5
6class scantable(sdtable):
7 """
8 The ASAP container for scans
9 """
[113]10
[484]11 def __init__(self, filename, average=None, unit=None):
[102]12 """
13 Create a scantable from a saved one or make a reference
14 Parameters:
[181]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
[102]22 scantable
[484]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.
[340]26 Over-rides the default selected by the reader
27 (input rpfits/sdfits/ms) or replaces the value
28 in existing scantables
[489]29 """
30 if average is None or type(average) is not bool:
31 average = rcParams['scantable.autoaverage']
32
[484]33 varlist = vars()
[226]34 self._vb = rcParams['verbose']
[113]35 self._p = None
[484]36
[181]37 if isinstance(filename,sdtable):
38 sdtable.__init__(self, filename)
[340]39 if unit is not None:
[346]40 self.set_fluxunit(unit)
[181]41 else:
[411]42 import os.path
43 if not os.path.exists(filename):
44 print "File '%s' not found." % (filename)
[226]45 return
[411]46 filename = os.path.expandvars(filename)
47 if os.path.isdir(filename):
[181]48 # crude check if asap table
[411]49 if os.path.exists(filename+'/table.info'):
[181]50 sdtable.__init__(self, filename)
[340]51 if unit is not None:
52 self.set_fluxunit(unit)
[181]53 else:
[411]54 print "The given file '%s'is not a valid asap table." % (filename)
[226]55 return
[484]56 else:
[181]57 from asap._asap import sdreader
[340]58 ifSel = -1
59 beamSel = -1
[345]60 r = sdreader(filename,ifSel,beamSel)
[181]61 print 'Importing data...'
[411]62 r._read([-1])
63 tbl = r._getdata()
[346]64 if unit is not None:
65 tbl.set_fluxunit(unit)
[489]66 if average:
67 from asap._asap import average as _av
[226]68 print 'Auto averaging integrations...'
[513]69 tbl2 = _av((tbl,),(),True,'none')
[226]70 sdtable.__init__(self,tbl2)
[513]71 del tbl2
[226]72 else:
73 sdtable.__init__(self,tbl)
[513]74 del r,tbl
[484]75 self._add_history("scantable", varlist)
[102]76
[451]77 def save(self, name=None, format=None, stokes=False, overwrite=False):
[116]78 """
[280]79 Store the scantable on disk. This can be an asap (aips++) Table, SDFITS,
80 Image FITS or MS2 format.
[116]81 Parameters:
[196]82 name: the name of the outputfile. For format="FITS" this
[489]83 is the directory file name into which all the files
[497]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)
[116]87 format: an optional file format. Default is ASAP.
[280]88 Allowed are - 'ASAP' (save as ASAP [aips++] Table),
[194]89 'SDFITS' (save as SDFITS file)
90 'FITS' (saves each row as a FITS Image)
[200]91 'ASCII' (saves as ascii text file)
[226]92 'MS2' (saves as an aips++
93 MeasurementSet V2)
[489]94 stokes: Convert to Stokes parameters (only available
95 currently with FITS and ASCII formats.
96 Default is False.
[411]97 overwrite: If the file should be overwritten if it exists.
[256]98 The default False is to return with warning
[411]99 without writing the output. USE WITH CARE.
[116]100 Example:
101 scan.save('myscan.asap')
102 scan.save('myscan.sdfits','SDFITS')
103 """
[411]104 from os import path
[226]105 if format is None: format = rcParams['scantable.save']
[256]106 suffix = '.'+format.lower()
107 if name is None or name =="":
108 name = 'scantable'+suffix
[283]109 print "No filename given. Using default name %s..." % name
[411]110 name = path.expandvars(name)
[256]111 if path.isfile(name) or path.isdir(name):
112 if not overwrite:
113 print "File %s already exists." % name
114 return
[451]115 format2 = format.upper()
116 if format2 == 'ASAP':
[116]117 self._save(name)
118 else:
119 from asap._asap import sdwriter as _sw
[451]120 w = _sw(format2)
[445]121 w.write(self, name, stokes)
[116]122 return
123
[102]124 def copy(self):
125 """
126 Return a copy of this scantable.
127 Parameters:
[113]128 none
[102]129 Example:
130 copiedscan = scan.copy()
131 """
[113]132 sd = scantable(sdtable._copy(self))
133 return sd
134
[102]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:
[513]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
[102]143 Example:
[513]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])
[102]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:
[113]155 s = sdtable._getsource(self,scanid)
[102]156 return scantable(s)
157 elif type(scanid) is int:
[381]158 s = sdtable._getscan(self,[scanid])
159 return scantable(s)
160 elif type(scanid) is list:
[113]161 s = sdtable._getscan(self,scanid)
[102]162 return scantable(s)
[381]163 else:
164 print "Illegal scanid type, use 'int' or 'list' if ints."
[102]165 except RuntimeError:
166 print "Couldn't find any match."
167
168 def __str__(self):
[381]169 return sdtable._summary(self,True)
[102]170
[381]171 def summary(self,filename=None, verbose=None):
[102]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
[381]177 verbose: print extra info such as the frequency table
178 The default (False) is taken from .asaprc
[102]179 """
[381]180 info = sdtable._summary(self, verbose)
181 if verbose is None: verbose = rcParams['scantable.verbosesummary']
[102]182 if filename is not None:
[256]183 if filename is "":
184 filename = 'scantable_summary.txt'
[415]185 from os.path import expandvars, isdir
[411]186 filename = expandvars(filename)
[415]187 if not isdir(filename):
[413]188 data = open(filename, 'w')
189 data.write(info)
190 data.close()
191 else:
192 print "Illegal file name '%s'." % (filename)
[102]193 print info
[484]194
[553]195 def set_cursor(self, beam=0, IF=0, pol=0):
[102]196 """
197 Set the spectrum for individual operations.
198 Parameters:
[553]199 beam, IF, pol: a number
[102]200 Example:
[256]201 scan.set_cursor(0,0,1)
[135]202 pol1sig = scan.stats(all=False) # returns std dev for beam=0
[181]203 # if=0, pol=1
[102]204 """
[484]205 varlist = vars()
[553]206 self.setbeam(beam)
207 self.setpol(pol)
208 self.setif(IF)
[484]209 self._add_history("set_cursor",varlist)
[102]210 return
211
[256]212 def get_cursor(self):
[113]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 """
[102]222 i = self.getbeam()
223 j = self.getif()
224 k = self.getpol()
[113]225 if self._vb:
[181]226 print "--------------------------------------------------"
[256]227 print " Cursor position"
[181]228 print "--------------------------------------------------"
[226]229 out = 'Beam=%d IF=%d Pol=%d ' % (i,j,k)
[113]230 print out
[102]231 return i,j,k
232
[436]233 def stats(self, stat='stddev', mask=None, allaxes=None):
[102]234 """
[135]235 Determine the specified statistic of the current beam/if/pol
[102]236 Takes a 'mask' as an optional parameter to specify which
237 channels should be excluded.
238 Parameters:
[135]239 stat: 'min', 'max', 'sumsq', 'sum', 'mean'
240 'var', 'stddev', 'avdev', 'rms', 'median'
241 mask: an optional mask specifying where the statistic
[102]242 should be determined.
[436]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)
[102]246 Example:
[113]247 scan.set_unit('channel')
[102]248 msk = scan.create_mask([100,200],[500,600])
[135]249 scan.stats(stat='mean', mask=m)
[102]250 """
[436]251 if allaxes is None: allaxes = rcParams['scantable.allaxes']
[135]252 from asap._asap import stats as _stats
[256]253 from numarray import array,zeros,Float
[102]254 if mask == None:
255 mask = ones(self.nchan())
[256]256 axes = ['Beam','IF','Pol','Time']
257
[436]258 beamSel,IFSel,polSel = (self.getbeam(),self.getif(),self.getpol())
259 if allaxes:
[256]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())]
[181]273 if self._vb:
[256]274 self._print_values(retval,stat,tm)
[436]275 self.setbeam(beamSel)
276 self.setif(IFSel)
277 self.setpol(polSel)
[256]278 return retval
[102]279
280 else:
[241]281 statval = _stats(self,mask,stat,-1)
[181]282 out = ''
283 for l in range(self.nrow()):
284 tm = self._gettime(l)
285 out += 'Time[%s]:\n' % (tm)
[436]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)
[181]289 out += '= %3.3f\n' % (statval[l])
290 out += "--------------------------------------------------\n"
[226]291
[181]292 if self._vb:
293 print "--------------------------------------------------"
294 print " ",stat
295 print "--------------------------------------------------"
296 print out
[256]297 retval = {'axes': axes, 'data': array(statval), 'cursor':(i,j,k)}
298 return retval
[102]299
[436]300 def stddev(self,mask=None, allaxes=None):
[135]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.
[436]308 allaxes: optional flag to show all or a cursor selected
[226]309 spectrum of Beam/IF/Pol. Default is all or taken
310 from .asaprc
[135]311
312 Example:
313 scan.set_unit('channel')
314 msk = scan.create_mask([100,200],[500,600])
315 scan.stddev(mask=m)
316 """
[436]317 if allaxes is None: allaxes = rcParams['scantable.allaxes']
318 return self.stats(stat='stddev',mask=mask, allaxes=allaxes);
[135]319
[436]320 def get_tsys(self, allaxes=None):
[113]321 """
322 Return the System temperatures.
323 Parameters:
[436]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)
[113]327 Returns:
328 a list of Tsys values.
329 """
[436]330 if allaxes is None: allaxes = rcParams['scantable.allaxes']
[256]331 from numarray import array,zeros,Float
332 axes = ['Beam','IF','Pol','Time']
333
[436]334 if allaxes:
[256]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())]
[181]348 if self._vb:
[256]349 self._print_values(retval,'Tsys',tm)
350 return retval
351
[102]352 else:
[256]353 i,j,k = (self.getbeam(),self.getif(),self.getpol())
354 statval = self._gettsys()
[181]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)
[256]362 out += '= %3.3f\n' % (statval[l])
363 out += "--------------------------------------------------\n"
364
[181]365 if self._vb:
366 print "--------------------------------------------------"
[256]367 print " TSys"
[181]368 print "--------------------------------------------------"
369 print out
[256]370 retval = {'axes': axes, 'data': array(statval), 'cursor':(i,j,k)}
371 return retval
[113]372
[407]373 def get_time(self, row=-1):
[113]374 """
375 Get a list of time stamps for the observations.
[181]376 Return a string for each integration in the scantable.
[113]377 Parameters:
[407]378 row: row no of integration. Default -1 return all rows
[113]379 Example:
380 none
381 """
382 out = []
[407]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)
[102]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'
[113]396 one of '*Hz','km/s','channel', ''
[102]397 """
[484]398 varlist = vars()
[113]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()
[484]405 self._add_history("set_unit",varlist)
[113]406
[484]407 def set_instrument(self, instr):
[358]408 """
409 Set the instrument for subsequent processing
410 Parameters:
[407]411 instr: Select from 'ATPKSMB', 'ATPKSHOH', 'ATMOPRA',
412 'DSS-43' (Tid), 'CEDUNA', and 'HOBART'
[358]413 """
414 self._setInstrument(instr)
[484]415 self._add_history("set_instument",vars())
[358]416
[276]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 """
[484]423 varlist = vars()
[276]424 inf = list(self._getcoordinfo())
425 inf[2] = doppler
426 self._setcoordinfo(inf)
427 if self._p: self.plot()
[484]428 self._add_history("set_doppler",vars())
[432]429
[226]430 def set_freqframe(self, frame=None):
[113]431 """
432 Set the frame type of the Spectral Axis.
433 Parameters:
[591]434 frame: an optional frame type, default 'LSRK'. Valid frames are:
435 'REST','TOPO','LSRD','LSRK','BARY',
436 'GEO','GALACTO','LGROUP','CMB'
[113]437 Examples:
438 scan.set_freqframe('BARY')
439 """
[484]440 if frame is None: frame = rcParams['scantable.freqframe']
441 varlist = vars()
[113]442 valid = ['REST','TOPO','LSRD','LSRK','BARY', \
443 'GEO','GALACTO','LGROUP','CMB']
[591]444
[113]445 if frame in valid:
446 inf = list(self._getcoordinfo())
447 inf[1] = frame
448 self._setcoordinfo(inf)
[484]449 self._add_history("set_freqframe",varlist)
[102]450 else:
[113]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
[102]464
[158]465 def get_abcissa(self, rowno=0):
[102]466 """
[158]467 Get the abcissa in the current coordinate setup for the currently
[113]468 selected Beam/IF/Pol
469 Parameters:
[226]470 rowno: an optional row number in the scantable. Default is the
471 first row, i.e. rowno=0
[113]472 Returns:
[407]473 The abcissa values and it's format string (as a dictionary)
[113]474 """
[256]475 abc = self._getabcissa(rowno)
[407]476 lbl = self._getabcissalabel(rowno)
[158]477 return abc, lbl
[407]478 #return {'abcissa':abc,'label':lbl}
[113]479
480 def create_mask(self, *args, **kwargs):
481 """
[102]482 Compute and return a mask based on [min,max] windows.
[189]483 The specified windows are to be INCLUDED, when the mask is
[113]484 applied.
[102]485 Parameters:
486 [min,max],[min2,max2],...
487 Pairs of start/end points specifying the regions
488 to be masked
[189]489 invert: optional argument. If specified as True,
490 return an inverted mask, i.e. the regions
491 specified are EXCLUDED
[513]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.
[102]495 Example:
[113]496 scan.set_unit('channel')
497
498 a)
[102]499 msk = scan.set_mask([400,500],[800,900])
[189]500 # masks everything outside 400 and 500
[113]501 # and 800 and 900 in the unit 'channel'
502
503 b)
504 msk = scan.set_mask([400,500],[800,900], invert=True)
[189]505 # masks the regions between 400 and 500
[113]506 # and 800 and 900 in the unit 'channel'
507
[102]508 """
[513]509 row = 0
510 if kwargs.has_key("row"):
511 row = kwargs.get("row")
512 data = self._getabcissa(row)
[113]513 u = self._getcoordinfo()[0]
514 if self._vb:
515 if u == "": u = "channel"
516 print "The current mask window unit is", u
[102]517 n = self.nchan()
[189]518 msk = zeros(n)
[484]519 for window in args:
[102]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):
[113]524 if data[i] >= window[0] and data[i] < window[1]:
[189]525 msk[i] = 1
[113]526 if kwargs.has_key('invert'):
527 if kwargs.get('invert'):
528 from numarray import logical_not
[513]529 msk = logical_not(msk)
[102]530 return msk
[113]531
[256]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())
[102]542
[402]543 def lines(self):
[391]544 """
[402]545 Print the list of known spectral lines
546 """
547 sdtable._lines(self)
548
[484]549 def set_restfreqs(self, freqs=None, unit='Hz', lines=None, source=None,
550 theif=None):
[402]551 """
[393]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.
[402]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()
[391]570 Parameters:
[402]571 freqs: list of rest frequencies
[393]572 unit: unit for rest frequency (default 'Hz')
[402]573 lines: list of known spectral lines names (alternative to freqs).
574 See possible list via scantable.lines()
[391]575 source: Source name (blank means all)
576 theif: IF (-1 means all)
577 Example:
[393]578 scan.set_restfreqs(freqs=1.4e9, source='NGC253', theif=2)
579 scan.set_restfreqs(freqs=[1.4e9,1.67e9])
[391]580 """
[484]581 varlist = vars()
[391]582 if source is None:
583 source = ""
584 if theif is None:
585 theif = -1
[393]586 t = type(freqs)
587 if t is int or t is float:
588 freqs = [freqs]
[402]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)
[484]597 self._add_history("set_restfreqs", varlist)
598
[391]599
600
[102]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 """
[407]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)
[484]620 self._add_history("flag_spectrum", vars())
[102]621 else:
622 print "Please specify a valid (Beam/IF/Pol)"
623 return
[113]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'
[135]631 col: which out of Beams/IFs/Pols should be colour stacked
[113]632 panel: set up multiple panels, currently not working.
633 """
[283]634 print "Warning! Not fully functional. Use plotter.plot() instead"
635
[113]636 validcol = {'Beam':self.nbeam(),'IF':self.nif(),'Pol':self.npol()}
[124]637
[113]638 validyax = ['spectrum','tsys']
[189]639 from asap.asaplot import ASAPlot
[113]640 if not self._p:
641 self._p = ASAPlot()
[189]642 #print "Plotting not enabled"
643 #return
644 if self._p.is_dead:
645 del self._p
646 self._p = ASAPlot()
[113]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
[124]654 self._p.hold()
655 self._p.clear()
[113]656 if panel == 'Time':
657 npan = self.nrow()
[124]658 self._p.set_panels(rows=npan)
[113]659 xlab,ylab,tlab = None,None,None
[226]660 self._vb = False
[256]661 sel = self.get_cursor()
[113]662 for i in range(npan):
[226]663 if npan > 1:
664 self._p.subplot(i)
[113]665 for j in range(validcol[col]):
666 x = None
667 y = None
668 m = None
669 tlab = self._getsourcename(i)
[124]670 import re
671 tlab = re.sub('_S','',tlab)
[113]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:
[158]687 x,xlab = self.get_abcissa(i)
[256]688 y = self._getspectrum(i)
[113]689 ylab = r'Flux'
[256]690 m = self._getmask(i)
[113]691 llab = col+' '+str(j)
692 self._p.set_line(label=llab)
693 self._p.plot(x,y,m)
[124]694 self._p.set_axes('xlabel',xlab)
695 self._p.set_axes('ylabel',ylab)
696 self._p.set_axes('title',tlab)
[113]697 self._p.release()
[256]698 self.set_cursor(sel[0],sel[1],sel[2])
[226]699 self._vb = rcParams['verbose']
[113]700 return
[256]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])
[484]719 out += "-"*80
720 out += "\n"
721 print "-"*80
[256]722 print " ", label
[484]723 print "-"*80
[256]724 print out
[484]725
726 def history(self):
727 hist = list(self._gethistory())
728 print "-"*80
729 for h in hist:
[489]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
[484]743 return
744
[513]745 #
746 # Maths business
747 #
748
[558]749 def average_time(self, mask=None, scanav=False, weight='tint'):
[513]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)
[558]759 scanav: True averages each scan separately
760 False (default) averages all scans together,
[524]761 weight: Weighting scheme. 'none', 'var' (1/var(spec)
[535]762 weighted), 'tsys' (1/Tsys**2 weighted), 'tint'
763 (integration time weighted) or 'tintsys' (Tint/Tsys**2).
764 The default is 'tint'
[513]765 Example:
766 # time average the scantable without using a mask
767 newscan = scan.average_time()
768 """
769 varlist = vars()
[524]770 if weight is None: weight = 'tint'
[513]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.
[532]992 weight: Weighting scheme. 'none' (default), 'var' (1/var(spec)
993 weighted), or 'tsys' (1/Tsys**2 weighted)
[513]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)
[579]1063 """
[513]1064 if insitu is None: insitu = rcParams['insitu']
1065 varlist = vars()
1066 if mask is None:
1067 from numarray import ones
[626]1068 mask = list(ones(self.nchan()))
[513]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:
[516]1106 scan2=scan.auto_poly_baseline(order=7)
[513]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):
[516]1115 raise RuntimeError, "Parameter 'edge' has to be an integer or a \
1116 pair of integers specified as a tuple"
[513]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 as _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
1258 def quotient(self, other, isreference=True, preserve=True):
1259 """
[582]1260 Return the quotient of a 'source' (on) scan and a 'reference' (off)
1261 scan.
[513]1262 The reference can have just one row, even if the signal has many.
1263 Otherwise they must have the same number of rows.
1264 The cursor of the output scan is set to 0
1265 Parameters:
[582]1266 other: the 'other' scan
[513]1267 isreference: if the 'other' scan is the reference (default)
1268 or source
1269 preserve: you can preserve (default) the continuum or
1270 remove it. The equations used are
[582]1271 preserve: Output = Toff * (on/off) - Toff
1272 remove: Output = Tref * (on/off) - Ton
[513]1273 Example:
1274 # src is a scantable for an 'on' scan, ref for an 'off' scantable
1275 q1 = src.quotient(ref)
1276 q2 = ref.quotient(src, isreference=False)
1277 # gives the same result
1278 """
1279 from asap._asap import quotient as _quot
1280 if isreference:
1281 return scantable(_quot(self, other, preserve))
1282 else:
1283 return scantable(_quot(other, self, preserve))
1284
1285 def __add__(self, other):
1286 varlist = vars()
1287 s = None
1288 if isinstance(other, scantable):
1289 from asap._asap import b_operate as _bop
1290 s = scantable(_bop(self, other, 'add', True))
1291 elif isinstance(other, float):
1292 from asap._asap import add as _add
1293 s = scantable(_add(self, other, True))
1294 else:
1295 print "Other input is not a scantable or float value"
1296 return
1297 s._add_history("operator +", varlist)
1298 return s
1299
1300
1301 def __sub__(self, other):
1302 """
1303 implicit on all axes and on Tsys
1304 """
1305 varlist = vars()
1306 s = None
1307 if isinstance(other, scantable):
1308 from asap._asap import b_operate as _bop
1309 s = scantable(_bop(self, other, 'sub', True))
1310 elif isinstance(other, float):
1311 from asap._asap import add as _add
1312 s = scantable(_add(self, -other, True))
1313 else:
1314 print "Other input is not a scantable or float value"
1315 return
1316 s._add_history("operator -", varlist)
1317 return s
1318
1319 def __mul__(self, other):
1320 """
1321 implicit on all axes and on Tsys
1322 """
1323 varlist = vars()
1324 s = None
1325 if isinstance(other, scantable):
1326 from asap._asap import b_operate as _bop
1327 s = scantable(_bop(self, other, 'mul', True))
1328 elif isinstance(other, float):
1329 if other == 0.0:
1330 print "Multiplying by zero is not recommended."
1331 return
1332 from asap._asap import scale as _sca
1333 s = scantable(_sca(self, other, True))
1334 else:
1335 print "Other input is not a scantable or float value"
1336 return
1337 s._add_history("operator *", varlist)
1338 return s
1339
1340
1341 def __div__(self, other):
1342 """
1343 implicit on all axes and on Tsys
1344 """
1345 varlist = vars()
1346 s = None
1347 if isinstance(other, scantable):
1348 from asap._asap import b_operate as _bop
1349 s = scantable(_bop(self, other, 'div', True))
1350 elif isinstance(other, float):
1351 if other == 0.0:
1352 print "Dividing by zero is not recommended"
1353 return
1354 from asap._asap import scale as _sca
1355 s = scantable(_sca(self, 1.0/other, True))
1356 else:
1357 print "Other input is not a scantable or float value"
1358 return
1359 s._add_history("operator /", varlist)
1360 return s
1361
[530]1362 def get_fit(self, row=0):
1363 """
1364 Print or return the stored fits for a row in the scantable
1365 Parameters:
1366 row: the row which the fit has been applied to.
1367 """
1368 if row > self.nrow():
1369 return
1370 from asap import asapfit
1371 fit = asapfit(self._getfit(row))
1372 if self._vb:
1373 print fit
1374 return
1375 else:
1376 return fit.as_dict()
1377
[484]1378 def _add_history(self, funcname, parameters):
1379 # create date
1380 sep = "##"
1381 from datetime import datetime
1382 dstr = datetime.now().strftime('%Y/%m/%d %H:%M:%S')
1383 hist = dstr+sep
1384 hist += funcname+sep#cdate+sep
1385 if parameters.has_key('self'): del parameters['self']
1386 for k,v in parameters.iteritems():
1387 if type(v) is dict:
1388 for k2,v2 in v.iteritems():
1389 hist += k2
1390 hist += "="
1391 if isinstance(v2,scantable):
1392 hist += 'scantable'
1393 elif k2 == 'mask':
[513]1394 if isinstance(v2,list) or isinstance(v2,tuple):
1395 hist += str(self._zip_mask(v2))
1396 else:
1397 hist += str(v2)
[484]1398 else:
[513]1399 hist += str(v2)
[484]1400 else:
1401 hist += k
1402 hist += "="
1403 if isinstance(v,scantable):
1404 hist += 'scantable'
1405 elif k == 'mask':
[513]1406 if isinstance(v,list) or isinstance(v,tuple):
1407 hist += str(self._zip_mask(v))
1408 else:
1409 hist += str(v)
[484]1410 else:
1411 hist += str(v)
1412 hist += sep
1413 hist = hist[:-2] # remove trailing '##'
1414 self._addhistory(hist)
1415
1416
1417 def _zip_mask(self, mask):
1418 mask = list(mask)
1419 i = 0
1420 segments = []
1421 while mask[i:].count(1):
1422 i += mask[i:].index(1)
1423 if mask[i:].count(0):
1424 j = i + mask[i:].index(0)
1425 else:
1426 j = len(mask)
1427 segments.append([i,j])
[513]1428 i = j
[484]1429 return segments
[626]1430 def _get_ordinate_label(self):
1431 fu = "("+self.get_fluxunit()+")"
1432 import re
1433 lbl = "Intensity"
1434 if re.match(".K.",fu):
1435 lbl = "Brightness Temperature "+ fu
1436 elif re.match(".Jy.",fu):
1437 lbl = "Flux density "+ fu
1438 return lbl
1439
Note: See TracBrowser for help on using the repository browser.