source: tags/asap0004/python/scantable.py@ 804

Last change on this file since 804 was 589, checked in by mar637, 19 years ago

Fix for asap0004

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 55.7 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 scan: The scantable
990 mask: An optional mask defining the region, where the
991 averaging will be applied. The output will have all
992 specified points masked.
993 weight: Weighting scheme. 'none' (default), 'var' (1/var(spec)
994 weighted), or 'tsys' (1/Tsys**2 weighted)
995 insitu: if False a new scantable is returned.
996 Otherwise, the scaling is done in-situ
997 The default is taken from .asaprc (False)
998 Example:
999 polav = average_pols(myscan)
1000 """
1001 if insitu is None: insitu = rcParams['insitu']
1002 varlist = vars()
1003
1004 if mask is None:
1005 mask = ()
1006 if not insitu:
1007 from asap._asap import averagepol as _avpol
1008 s = scantable(_avpol(self, mask, weight))
1009 s._add_history("average_pol",varlist)
1010 return s
1011 else:
1012 from asap._asap import averagepol_insitu as _avpol
1013 _avpol(self, mask, weight)
1014 self._add_history("average_pol",varlist)
1015 return
1016
1017 def smooth(self, kernel="hanning", width=5.0, insitu=None, allaxes=None):
1018 """
1019 Smooth the spectrum by the specified kernel (conserving flux).
1020 Parameters:
1021 scan: The input scan
1022 kernel: The type of smoothing kernel. Select from
1023 'hanning' (default), 'gaussian' and 'boxcar'.
1024 The first three characters are sufficient.
1025 width: The width of the kernel in pixels. For hanning this is
1026 ignored otherwise it defauls to 5 pixels.
1027 For 'gaussian' it is the Full Width Half
1028 Maximum. For 'boxcar' it is the full width.
1029 insitu: if False a new scantable is returned.
1030 Otherwise, the scaling is done in-situ
1031 The default is taken from .asaprc (False)
1032 allaxes: If True (default) apply to all spectra. Otherwise
1033 apply only to the selected (beam/pol/if)spectra only
1034 The default is taken from .asaprc (True if none)
1035 Example:
1036 none
1037 """
1038 if allaxes is None: allaxes = rcParams['scantable.allaxes']
1039 if insitu is None: insitu = rcParams['insitu']
1040 varlist = vars()
1041 if not insitu:
1042 from asap._asap import smooth as _smooth
1043 s = scantable(_smooth(self,kernel,width,allaxes))
1044 s._add_history("smooth", varlist)
1045 return s
1046 else:
1047 from asap._asap import smooth_insitu as _smooth
1048 _smooth(self,kernel,width,allaxes)
1049 self._add_history("smooth", varlist)
1050 return
1051
1052 def poly_baseline(self, mask=None, order=0, insitu=None):
1053 """
1054 Return a scan which has been baselined (all rows) by a polynomial.
1055 Parameters:
1056 scan: a scantable
1057 mask: an optional mask
1058 order: the order of the polynomial (default is 0)
1059 insitu: if False a new scantable is returned.
1060 Otherwise, the scaling is done in-situ
1061 The default is taken from .asaprc (False)
1062 Example:
1063 # return a scan baselined by a third order polynomial,
1064 # not using a mask
1065 bscan = scan.poly_baseline(order=3)
1066 """
1067 if insitu is None: insitu = rcParams['insitu']
1068 varlist = vars()
1069 if mask is None:
1070 from numarray import ones
1071 mask = list(ones(scan.nchan()))
1072 from asap.asapfitter import fitter
1073 f = fitter()
1074 f._verbose(True)
1075 f.set_scan(self, mask)
1076 f.set_function(poly=order)
1077 sf = f.auto_fit(insitu)
1078 if insitu:
1079 self._add_history("poly_baseline", varlist)
1080 return
1081 else:
1082 sf._add_history("poly_baseline", varlist)
1083 return sf
1084
1085 def auto_poly_baseline(self, mask=None, edge=(0,0), order=0,
1086 threshold=3,insitu=None):
1087 """
1088 Return a scan which has been baselined (all rows) by a polynomial.
1089 Spectral lines are detected first using linefinder and masked out
1090 to avoid them affecting the baseline solution.
1091
1092 Parameters:
1093 scan: a scantable
1094 mask: an optional mask retreived from scantable
1095 edge: an optional number of channel to drop at
1096 the edge of spectrum. If only one value is
1097 specified, the same number will be dropped from
1098 both sides of the spectrum. Default is to keep
1099 all channels
1100 order: the order of the polynomial (default is 0)
1101 threshold: the threshold used by line finder. It is better to
1102 keep it large as only strong lines affect the
1103 baseline solution.
1104 insitu: if False a new scantable is returned.
1105 Otherwise, the scaling is done in-situ
1106 The default is taken from .asaprc (False)
1107
1108 Example:
1109 scan2=scan.auto_poly_baseline(order=7)
1110 """
1111 if insitu is None: insitu = rcParams['insitu']
1112 varlist = vars()
1113 from asap.asapfitter import fitter
1114 from asap.asaplinefind import linefinder
1115 from asap import _is_sequence_or_number as _is_valid
1116
1117 if not _is_valid(edge, int):
1118 raise RuntimeError, "Parameter 'edge' has to be an integer or a \
1119 pair of integers specified as a tuple"
1120
1121 # setup fitter
1122 f = fitter()
1123 f._verbose(True)
1124 f.set_function(poly=order)
1125
1126 # setup line finder
1127 fl=linefinder()
1128 fl.set_options(threshold=threshold)
1129
1130 if not insitu:
1131 workscan=self.copy()
1132 else:
1133 workscan=self
1134
1135 vb=workscan._vb
1136 # remember the verbose parameter and selection
1137 workscan._vb=False
1138 sel=workscan.get_cursor()
1139 rows=range(workscan.nrow())
1140 for i in range(workscan.nbeam()):
1141 workscan.setbeam(i)
1142 for j in range(workscan.nif()):
1143 workscan.setif(j)
1144 for k in range(workscan.npol()):
1145 workscan.setpol(k)
1146 if f._vb:
1147 print "Processing:"
1148 print 'Beam[%d], IF[%d], Pol[%d]' % (i,j,k)
1149 for iRow in rows:
1150 fl.set_scan(workscan,mask,edge)
1151 fl.find_lines(iRow)
1152 f.set_scan(workscan, fl.get_mask())
1153 f.x=workscan._getabcissa(iRow)
1154 f.y=workscan._getspectrum(iRow)
1155 f.data=None
1156 f.fit()
1157 x=f.get_parameters()
1158 workscan._setspectrum(f.fitter.getresidual(),iRow)
1159 workscan.set_cursor(sel[0],sel[1],sel[2])
1160 workscan._vb = vb
1161 if not insitu:
1162 return workscan
1163
1164 def rotate_linpolphase(self, angle, allaxes=None):
1165 """
1166 Rotate the phase of the complex polarization O=Q+iU correlation.
1167 This is always done in situ in the raw data. So if you call this
1168 function more than once then each call rotates the phase further.
1169 Parameters:
1170 angle: The angle (degrees) to rotate (add) by.
1171 allaxes: If True apply to all spectra. Otherwise
1172 apply only to the selected (beam/pol/if)spectra only.
1173 The default is taken from .asaprc (True if none)
1174 Examples:
1175 scan.rotate_linpolphase(2.3)
1176 """
1177 if allaxes is None: allaxes = rcParams['scantable.allaxes']
1178 varlist = vars()
1179 from asap._asap import _rotate_linpolphase as _rotate
1180 _rotate(self, angle, allaxes)
1181 self._add_history("rotate_linpolphase", varlist)
1182 return
1183
1184
1185 def rotate_xyphase(self, angle, allaxes=None):
1186 """
1187 Rotate the phase of the XY correlation. This is always done in situ
1188 in the data. So if you call this function more than once
1189 then each call rotates the phase further.
1190 Parameters:
1191 angle: The angle (degrees) to rotate (add) by.
1192 allaxes: If True apply to all spectra. Otherwise
1193 apply only to the selected (beam/pol/if)spectra only.
1194 The default is taken from .asaprc (True if none)
1195 Examples:
1196 scan.rotate_xyphase(2.3)
1197 """
1198 if allaxes is None: allaxes = rcParams['scantable.allaxes']
1199 varlist = vars()
1200 from asap._asap import rotate_xyphase as _rotate_xyphase
1201 _rotate_xyphase(self, angle, allaxes)
1202 self._add_history("rotate_xyphase", varlist)
1203 return
1204
1205
1206 def add(self, offset, insitu=None, allaxes=None):
1207 """
1208 Return a scan where all spectra have the offset added
1209 Parameters:
1210 offset: the offset
1211 insitu: if False a new scantable is returned.
1212 Otherwise, the scaling is done in-situ
1213 The default is taken from .asaprc (False)
1214 allaxes: if True apply to all spectra. Otherwise
1215 apply only to the selected (beam/pol/if)spectra only
1216 The default is taken from .asaprc (True if none)
1217 """
1218 if allaxes is None: allaxes = rcParams['scantable.allaxes']
1219 if insitu is None: insitu = rcParams['insitu']
1220 varlist = vars()
1221 if not insitu:
1222 from asap._asap import add as _add
1223 s = scantable(_add(self, offset, allaxes))
1224 s._add_history("add",varlist)
1225 return s
1226 else:
1227 from asap._asap import add_insitu as _add
1228 _add(self, offset, allaxes)
1229 self._add_history("add",varlist)
1230 return
1231
1232 def scale(self, factor, insitu=None, allaxes=None, tsys=True):
1233 """
1234 Return a scan where all spectra are scaled by the give 'factor'
1235 Parameters:
1236 factor: the scaling factor
1237 insitu: if False a new scantable is returned.
1238 Otherwise, the scaling is done in-situ
1239 The default is taken from .asaprc (False)
1240 allaxes: if True apply to all spectra. Otherwise
1241 apply only to the selected (beam/pol/if)spectra only.
1242 The default is taken from .asaprc (True if none)
1243 tsys: if True (default) then apply the operation to Tsys
1244 as well as the data
1245 """
1246 if allaxes is None: allaxes = rcParams['scantable.allaxes']
1247 if insitu is None: insitu = rcParams['insitu']
1248 varlist = vars()
1249 if not insitu:
1250 from asap._asap import scale as _scale
1251 s = scantable(_scale(self, factor, allaxes, tsys))
1252 s._add_history("scale",varlist)
1253 return s
1254 else:
1255 from asap._asap import scale_insitu as _scale
1256 _scale(self, factor, allaxes, tsys)
1257 self._add_history("scale",varlist)
1258 return
1259
1260
1261 def quotient(self, other, isreference=True, preserve=True):
1262 """
1263 Return the quotient of a 'source' (on) scan and a 'reference' (off)
1264 scan.
1265 The reference can have just one row, even if the signal has many.
1266 Otherwise they must have the same number of rows.
1267 The cursor of the output scan is set to 0
1268 Parameters:
1269 other: the 'other' scan
1270 isreference: if the 'other' scan is the reference (default)
1271 or source
1272 preserve: you can preserve (default) the continuum or
1273 remove it. The equations used are
1274 preserve: Output = Toff * (on/off) - Toff
1275 remove: Output = Tref * (on/off) - Ton
1276 Example:
1277 # src is a scantable for an 'on' scan, ref for an 'off' scantable
1278 q1 = src.quotient(ref)
1279 q2 = ref.quotient(src, isreference=False)
1280 # gives the same result
1281 """
1282 from asap._asap import quotient as _quot
1283 if isreference:
1284 return scantable(_quot(self, other, preserve))
1285 else:
1286 return scantable(_quot(other, self, preserve))
1287
1288 def __add__(self, other):
1289 varlist = vars()
1290 s = None
1291 if isinstance(other, scantable):
1292 from asap._asap import b_operate as _bop
1293 s = scantable(_bop(self, other, 'add', True))
1294 elif isinstance(other, float):
1295 from asap._asap import add as _add
1296 s = scantable(_add(self, other, True))
1297 else:
1298 print "Other input is not a scantable or float value"
1299 return
1300 s._add_history("operator +", varlist)
1301 return s
1302
1303
1304 def __sub__(self, other):
1305 """
1306 implicit on all axes and on Tsys
1307 """
1308 varlist = vars()
1309 s = None
1310 if isinstance(other, scantable):
1311 from asap._asap import b_operate as _bop
1312 s = scantable(_bop(self, other, 'sub', True))
1313 elif isinstance(other, float):
1314 from asap._asap import add as _add
1315 s = scantable(_add(self, -other, True))
1316 else:
1317 print "Other input is not a scantable or float value"
1318 return
1319 s._add_history("operator -", varlist)
1320 return s
1321
1322 def __mul__(self, other):
1323 """
1324 implicit on all axes and on Tsys
1325 """
1326 varlist = vars()
1327 s = None
1328 if isinstance(other, scantable):
1329 from asap._asap import b_operate as _bop
1330 s = scantable(_bop(self, other, 'mul', True))
1331 elif isinstance(other, float):
1332 if other == 0.0:
1333 print "Multiplying by zero is not recommended."
1334 return
1335 from asap._asap import scale as _sca
1336 s = scantable(_sca(self, other, True))
1337 else:
1338 print "Other input is not a scantable or float value"
1339 return
1340 s._add_history("operator *", varlist)
1341 return s
1342
1343
1344 def __div__(self, other):
1345 """
1346 implicit on all axes and on Tsys
1347 """
1348 varlist = vars()
1349 s = None
1350 if isinstance(other, scantable):
1351 from asap._asap import b_operate as _bop
1352 s = scantable(_bop(self, other, 'div', True))
1353 elif isinstance(other, float):
1354 if other == 0.0:
1355 print "Dividing by zero is not recommended"
1356 return
1357 from asap._asap import scale as _sca
1358 s = scantable(_sca(self, 1.0/other, True))
1359 else:
1360 print "Other input is not a scantable or float value"
1361 return
1362 s._add_history("operator /", varlist)
1363 return s
1364
1365 def get_fit(self, row=0):
1366 """
1367 Print or return the stored fits for a row in the scantable
1368 Parameters:
1369 row: the row which the fit has been applied to.
1370 """
1371 if row > self.nrow():
1372 return
1373 from asap import asapfit
1374 fit = asapfit(self._getfit(row))
1375 if self._vb:
1376 print fit
1377 return
1378 else:
1379 return fit.as_dict()
1380
1381 def _add_history(self, funcname, parameters):
1382 # create date
1383 sep = "##"
1384 from datetime import datetime
1385 dstr = datetime.now().strftime('%Y/%m/%d %H:%M:%S')
1386 hist = dstr+sep
1387 hist += funcname+sep#cdate+sep
1388 if parameters.has_key('self'): del parameters['self']
1389 for k,v in parameters.iteritems():
1390 if type(v) is dict:
1391 for k2,v2 in v.iteritems():
1392 hist += k2
1393 hist += "="
1394 if isinstance(v2,scantable):
1395 hist += 'scantable'
1396 elif k2 == 'mask':
1397 if isinstance(v2,list) or isinstance(v2,tuple):
1398 hist += str(self._zip_mask(v2))
1399 else:
1400 hist += str(v2)
1401 else:
1402 hist += str(v2)
1403 else:
1404 hist += k
1405 hist += "="
1406 if isinstance(v,scantable):
1407 hist += 'scantable'
1408 elif k == 'mask':
1409 if isinstance(v,list) or isinstance(v,tuple):
1410 hist += str(self._zip_mask(v))
1411 else:
1412 hist += str(v)
1413 else:
1414 hist += str(v)
1415 hist += sep
1416 hist = hist[:-2] # remove trailing '##'
1417 self._addhistory(hist)
1418
1419
1420 def _zip_mask(self, mask):
1421 mask = list(mask)
1422 i = 0
1423 segments = []
1424 while mask[i:].count(1):
1425 i += mask[i:].index(1)
1426 if mask[i:].count(0):
1427 j = i + mask[i:].index(0)
1428 else:
1429 j = len(mask)
1430 segments.append([i,j])
1431 i = j
1432 return segments
Note: See TracBrowser for help on using the repository browser.