source: trunk/python/asapmath.py@ 502

Last change on this file since 502 was 489, checked in by mar637, 20 years ago
  • added history to all functions, which modify the data
  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 18.5 KB
Line 
1from scantable import scantable
2from asap import rcParams
3
4def average_time(*args, **kwargs):
5 """
6 Return the (time) average of a scan or list of scans. [in channels only]
7 The cursor of the output scan is set to 0
8 Parameters:
9 one scan or comma separated scans
10 mask: an optional mask (only used for 'var' and 'tsys' weighting)
11 scanav: True (default) averages each scan separately
12 False averages all scans together,
13 weight: Weighting scheme. 'none' (default), 'var' (variance
14 weighted), 'tsys'
15 Example:
16 # return a time averaged scan from scana and scanb
17 # without using a mask
18 scanav = average_time(scana,scanb)
19 # return the (time) averaged scan, i.e. the average of
20 # all correlator cycles
21 scanav = average_time(scan)
22
23 """
24 scanAv = True
25 if kwargs.has_key('scanav'):
26 scanAv = kwargs.get('scanav')
27 weight = 'none'
28 if kwargs.has_key('weight'):
29 weight = kwargs.get('weight')
30 mask = ()
31 if kwargs.has_key('mask'):
32 mask = kwargs.get('mask')
33 varlist = vars()
34 lst = tuple(args)
35 del varlist["kwargs"]
36 varlist["args"] = "%d scantables" % len(lst)
37 # need special formatting her for history...
38
39 from asap._asap import average as _av
40 for s in lst:
41 if not isinstance(s,scantable):
42 print "Please give a list of scantables"
43 return
44 s = scantable(_av(lst, mask, scanAv, weight))
45 s._add_history("average_time",varlist)
46 return s
47
48def quotient(source, reference, preserve=True):
49 """
50 Return the quotient of a 'source' (signal) scan and a 'reference' scan.
51 The reference can have just one row, even if the signal has many. Otherwise
52 they must have the same number of rows.
53 The cursor of the output scan is set to 0
54 Parameters:
55 source: the 'on' scan
56 reference: the 'off' scan
57 preserve: you can preserve (default) the continuum or
58 remove it. The equations used are
59 preserve - Output = Tref * (sig/ref) - Tref
60 remove - Output = Tref * (sig/ref) - Tsig
61 """
62 from asap._asap import quotient as _quot
63 return scantable(_quot(source, reference, preserve))
64
65def simple_math(left, right, op='add', tsys=True):
66 """
67 Apply simple mathematical binary operations to two
68 scan tables, returning the result in a new scan table.
69 The operation is applied to both the correlations and the TSys data
70 The cursor of the output scan is set to 0
71 Parameters:
72 left: the 'left' scan
73 right: the 'right' scan
74 op: the operation: 'add' (default), 'sub', 'mul', 'div'
75 tsys: if True (default) then apply the operation to Tsys
76 as well as the data
77 """
78 varlist = vars()
79 if not isinstance(left,scantable) and not isinstance(right,scantable):
80 print "Please provide two scantables as input"
81 return
82 from asap._asap import b_operate as _bop
83 s = scantable(_bop(left, right, op, tsys))
84 s._add_history("simple_math", varlist)
85 return s
86
87def scale(scan, factor, insitu=None, allaxes=None, tsys=True):
88 """
89 Return a scan where all spectra are scaled by the give 'factor'
90 Parameters:
91 scan: a scantable
92 factor: the scaling factor
93 insitu: if False a new scantable is returned.
94 Otherwise, the scaling is done in-situ
95 The default is taken from .asaprc (False)
96 allaxes: if True apply to all spectra. Otherwise
97 apply only to the selected (beam/pol/if)spectra only.
98 The default is taken from .asaprc (True if none)
99 tsys: if True (default) then apply the operation to Tsys
100 as well as the data
101 """
102 if allaxes is None: allaxes = rcParams['scantable.allaxes']
103 if insitu is None: insitu = rcParams['insitu']
104 varlist = vars()
105 if not insitu:
106 from asap._asap import scale as _scale
107 s = scantable(_scale(scan, factor, allaxes, tsys))
108 s._add_history("scale",varlist)
109 return s
110 else:
111 from asap._asap import scale_insitu as _scale
112 _scale(scan, factor, allaxes, tsys)
113 scan._add_history("scale",varlist)
114 return
115
116def add(scan, offset, insitu=None, allaxes=None):
117 """
118 Return a scan where all spectra have the offset added
119 Parameters:
120 scan: a scantable
121 offset: the offset
122 insitu: if False a new scantable is returned.
123 Otherwise, the scaling is done in-situ
124 The default is taken from .asaprc (False)
125 allaxes: if True apply to all spectra. Otherwise
126 apply only to the selected (beam/pol/if)spectra only
127 The default is taken from .asaprc (True if none)
128 """
129 if allaxes is None: allaxes = rcParams['scantable.allaxes']
130 if insitu is None: insitu = rcParams['insitu']
131 if not insitu:
132 from asap._asap import add as _add
133 return scantable(_add(scan, offset, allaxes))
134 else:
135 from asap._asap import add_insitu as _add
136 _add(scan, offset, allaxes)
137 return
138
139def convert_flux(scan, jyperk=None, eta=None, d=None, insitu=None,
140 allaxes=None):
141 """
142 Return a scan where all spectra are converted to either Jansky or Kelvin
143 depending upon the flux units of the scan table. By default the
144 function tries to look the values up internally. If it can't find
145 them (or if you want to over-ride), you must specify EITHER jyperk
146 OR eta (and D which it will try to look up also if you don't
147 set it). jyperk takes precedence if you set both.
148 Parameters:
149 scan: a scantable
150 jyperk: the Jy / K conversion factor
151 eta: the aperture efficiency
152 d: the geomtric diameter (metres)
153 insitu: if False a new scantable is returned.
154 Otherwise, the scaling is done in-situ
155 The default is taken from .asaprc (False)
156 allaxes: if True apply to all spectra. Otherwise
157 apply only to the selected (beam/pol/if)spectra only
158 The default is taken from .asaprc (True if none)
159 """
160 if allaxes is None: allaxes = rcParams['scantable.allaxes']
161 if insitu is None: insitu = rcParams['insitu']
162 varlist = vars()
163 if jyperk is None: jyperk = -1.0
164 if d is None: d = -1.0
165 if eta is None: eta = -1.0
166 if not insitu:
167 from asap._asap import convertflux as _convert
168 s = scantable(_convert(scan, d, eta, jyperk, allaxes))
169 s._add_history("convert_flux", varlist)
170 return s
171 else:
172 from asap._asap import convertflux_insitu as _convert
173 _convert(scan, d, eta, jyperk, allaxes)
174 scan._add_history("convert_flux", varlist)
175 return
176
177def gain_el(scan, poly=None, filename="", method="linear",
178 insitu=None, allaxes=None):
179 """
180 Return a scan after applying a gain-elevation correction. The correction
181 can be made via either a polynomial or a table-based interpolation
182 (and extrapolation if necessary).
183 You specify polynomial coefficients, an ascii table or neither.
184 If you specify neither, then a polynomial correction will be made
185 with built in coefficients known for certain telescopes (an error will
186 occur if the instrument is not known). The data and Tsys are *divided*
187 by the scaling factors.
188 Parameters:
189 scan: a scantable
190 poly: Polynomial coefficients (default None) to compute a
191 gain-elevation correction as a function of
192 elevation (in degrees).
193 filename: The name of an ascii file holding correction factors.
194 The first row of the ascii file must give the column
195 names and these MUST include columns
196 "ELEVATION" (degrees) and "FACTOR" (multiply data by this)
197 somewhere.
198 The second row must give the data type of the column. Use
199 'R' for Real and 'I' for Integer. An example file
200 would be (actual factors are arbitrary) :
201
202 TIME ELEVATION FACTOR
203 R R R
204 0.1 0 0.8
205 0.2 20 0.85
206 0.3 40 0.9
207 0.4 60 0.85
208 0.5 80 0.8
209 0.6 90 0.75
210 method: Interpolation method when correcting from a table. Values
211 are "nearest", "linear" (default), "cubic" and "spline"
212 insitu: if False a new scantable is returned.
213 Otherwise, the scaling is done in-situ
214 The default is taken from .asaprc (False)
215 allaxes: if True apply to all spectra. Otherwise
216 apply only to the selected (beam/pol/if) spectra only
217 The default is taken from .asaprc (True if none)
218 """
219 if allaxes is None: allaxes = rcParams['scantable.allaxes']
220 if insitu is None: insitu = rcParams['insitu']
221 varlist = vars()
222 if poly is None:
223 poly = ()
224 from os.path import expandvars
225 filename = expandvars(filename)
226 if not insitu:
227 from asap._asap import gainel as _gainEl
228 s = scantable(_gainEl(scan, poly, filename, method, allaxes))
229 s._add_history("gain_el", varlist)
230 return s
231 else:
232 from asap._asap import gainel_insitu as _gainEl
233 _gainEl(scan, poly, filename, method, allaxes)
234 scan._add_history("gain_el", varlist)
235 return
236
237def freq_align(scan, reftime=None, method='cubic', perif=False, insitu=None):
238 """
239 Return a scan where all rows have been aligned in frequency. The
240 alignment frequency frame (e.g. LSRK) is that set by function
241 set_freqframe.
242 scan: a scantable
243 reftime: reference time to align at. By default, the time of
244 the first row of data is used.
245 method: Interpolation method for regridding the spectra. Choose
246 from "nearest", "linear", "cubic" (default) and "spline"
247 perif: Generate aligners per freqID (no doppler tracking) or
248 per IF (scan-based doppler tracking)
249 insitu: if False a new scantable is returned.
250 Otherwise, the scaling is done in-situ
251 The default is taken from .asaprc (False)
252 """
253 if insitu is None: insitu = rcParams['insitu']
254 varlist = vars()
255 if reftime is None: reftime = ''
256 perfreqid = not perif
257 if not insitu:
258 from asap._asap import freq_align as _align
259 s = scantable(_align(scan, reftime, method, perfreqid))
260 s._add_history("freq_align", varlist)
261 return s
262 else:
263 from asap._asap import freq_align_insitu as _align
264 _align(scan, reftime, method, perfreqid)
265 scan._add_history("freq_align", varlist)
266 return
267
268def opacity(scan, tau, insitu=None, allaxes=None):
269 """
270 Return a scan after applying an opacity correction. The data
271 and Tsys are multiplied by the correction factor.
272 Parameters:
273 scan: a scantable
274 tau: Opacity from which the correction factor is exp(tau*ZD)
275 where ZD is the zenith-distance
276 insitu: if False a new scantable is returned.
277 Otherwise, the scaling is done in-situ
278 The default is taken from .asaprc (False)
279 allaxes: if True apply to all spectra. Otherwise
280 apply only to the selected (beam/pol/if)spectra only
281 The default is taken from .asaprc (True if none)
282 """
283 if allaxes is None: allaxes = rcParams['scantable.allaxes']
284 if insitu is None: insitu = rcParams['insitu']
285 varlist = vars()
286 if not insitu:
287 from asap._asap import opacity as _opacity
288 s = scantable(_opacity(scan, tau, allaxes))
289 s._add_history("opacity", varlist)
290 return s
291 else:
292 from asap._asap import opacity_insitu as _opacity
293 _opacity(scan, tau, allaxes)
294 scan._add_history("opacity", varlist)
295 return
296
297def bin(scan, width=5, insitu=None):
298 """
299 Return a scan where all spectra have been binned up
300 width: The bin width (default=5) in pixels
301 insitu: if False a new scantable is returned.
302 Otherwise, the scaling is done in-situ
303 The default is taken from .asaprc (False)
304 """
305 if insitu is None: insitu = rcParams['insitu']
306 varlist = vars()
307 if not insitu:
308 from asap._asap import bin as _bin
309 s = scantable(_bin(scan, width))
310 s._add_history("bin",varlist)
311 return s
312 else:
313 from asap._asap import bin_insitu as _bin
314 _bin(scan, width)
315 scan._add_history("bin",varlist)
316 return
317
318def resample(scan, width=5, method='cubic', insitu=None):
319 """
320 Return a scan where all spectra have been binned up
321 width: The bin width (default=5) in pixels
322 method: Interpolation method when correcting from a table. Values
323 are "nearest", "linear", "cubic" (default) and "spline"
324 insitu: if False a new scantable is returned.
325 Otherwise, the scaling is done in-situ
326 The default is taken from .asaprc (False)
327 """
328 if insitu is None: insitu = rcParams['insitu']
329 varlist = vars()
330 if not insitu:
331 from asap._asap import resample as _resample
332 s = scantable(_resample(scan, method, width))
333 s._add_history("resample",varlist)
334 return s
335 else:
336 from asap._asap import resample_insitu as _resample
337 _resample(scan, method, width)
338 scan._add_history("resample",varlist)
339 return
340
341def average_pol(scan, mask=None, weight='none', insitu=None):
342 """
343 Average the Polarisations together.
344 The polarisation cursor of the output scan is set to 0
345 Parameters:
346 scan: The scantable
347 mask: An optional mask defining the region, where the
348 averaging will be applied. The output will have all
349 specified points masked.
350 weight: Weighting scheme. 'none' (default), or 'var' (variance
351 weighted)
352 insitu: if False a new scantable is returned.
353 Otherwise, the scaling is done in-situ
354 The default is taken from .asaprc (False)
355 Example:
356 polav = average_pols(myscan)
357 """
358 if insitu is None: insitu = rcParams['insitu']
359 varlist = vars()
360
361 if mask is None:
362 mask = ()
363 if not insitu:
364 from asap._asap import averagepol as _avpol
365 s = scantable(_avpol(scan, mask, weight))
366 s._add_history("average_pol",varlist)
367 return s
368 else:
369 from asap._asap import averagepol_insitu as _avpol
370 _avpol(scan, mask, weight)
371 scan._add_history("average_pol",varlist)
372 return
373
374def smooth(scan, kernel="hanning", width=5.0, insitu=None, allaxes=None):
375 """
376 Smooth the spectrum by the specified kernel (conserving flux).
377 Parameters:
378 scan: The input scan
379 kernel: The type of smoothing kernel. Select from
380 'hanning' (default), 'gaussian' and 'boxcar'.
381 The first three characters are sufficient.
382 width: The width of the kernel in pixels. For hanning this is
383 ignored otherwise it defauls to 5 pixels.
384 For 'gaussian' it is the Full Width Half
385 Maximum. For 'boxcar' it is the full width.
386 insitu: if False a new scantable is returned.
387 Otherwise, the scaling is done in-situ
388 The default is taken from .asaprc (False)
389 allaxes: If True (default) apply to all spectra. Otherwise
390 apply only to the selected (beam/pol/if)spectra only
391 The default is taken from .asaprc (True if none)
392 Example:
393 none
394 """
395 if allaxes is None: allaxes = rcParams['scantable.allaxes']
396 if insitu is None: insitu = rcParams['insitu']
397 varlist = vars()
398 if not insitu:
399 from asap._asap import smooth as _smooth
400 s = scantable(_smooth(scan,kernel,width,allaxes))
401 s._add_history("smooth", varlist)
402 return s
403 else:
404 from asap._asap import smooth_insitu as _smooth
405 _smooth(scan,kernel,width,allaxes)
406 scan._add_history("smooth", varlist)
407 return
408
409def poly_baseline(scan, mask=None, order=0, insitu=None):
410 """
411 Return a scan which has been baselined (all rows) by a polynomial.
412 Parameters:
413 scan: a scantable
414 mask: an optional mask
415 order: the order of the polynomial (default is 0)
416 insitu: if False a new scantable is returned.
417 Otherwise, the scaling is done in-situ
418 The default is taken from .asaprc (False)
419 Example:
420 # return a scan baselined by a third order polynomial,
421 # not using a mask
422 bscan = poly_baseline(scan, order=3)
423 """
424 if insitu is None: insitu = rcParams['insitu']
425 varlist = vars()
426 from asap.asapfitter import fitter
427 if mask is None:
428 from numarray import ones
429 mask = tuple(ones(scan.nchan()))
430 f = fitter()
431 f._verbose(True)
432 f.set_scan(scan, mask)
433 f.set_function(poly=order)
434 sf = f.auto_fit(insitu)
435 sf._add_history("poly_baseline", varlist)
436 return sf
437
438def rotate_xyphase(scan, angle, allaxes=None):
439 """
440 Rotate the phase of the XY correlation. This is always done in situ
441 in the data. So if you call this function more than once
442 then each call rotates the phase further.
443 Parameters:
444 angle: The angle (degrees) to rotate (add) by.
445 allaxes: If True apply to all spectra. Otherwise
446 apply only to the selected (beam/pol/if)spectra only.
447 The default is taken from .asaprc (True if none)
448 Examples:
449 rotate_xyphase(scan, 2.3)
450 """
451 if allaxes is None: allaxes = rcParams['scantable.allaxes']
452 varlist = vars()
453 from asap._asap import rotate_xyphase as _rotate_xyphase
454 _rotate_xyphase(scan, angle, allaxes)
455 s._add_history("rotate_xyphase", varlist)
456 return
Note: See TracBrowser for help on using the repository browser.