source: branches/alma/python/asapmath.py@ 1620

Last change on this file since 1620 was 1614, checked in by Takeshi Nakazato, 15 years ago

New Development: No

JIRA Issue: Yes CAS-729, CAS-1147

Ready to Release: Yes

Interface Changes: No

What Interface Changed: Please list interface changes

Test Programs: List test programs

Put in Release Notes: No

Module(s): Module Names change impacts.

Description: Describe your changes here...

  1. Added level parameter to print_log()
  2. Replaced casalog.post() to asaplog.push() + print_log().


  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 14.5 KB
Line 
1from asap.scantable import scantable
2from asap import rcParams
3from asap import print_log
4from asap import selector
5from asap import asaplog
6
7def average_time(*args, **kwargs):
8 """
9 Return the (time) average of a scan or list of scans. [in channels only]
10 The cursor of the output scan is set to 0
11 Parameters:
12 one scan or comma separated scans or a list of scans
13 mask: an optional mask (only used for 'var' and 'tsys' weighting)
14 scanav: True averages each scan separately.
15 False (default) averages all scans together,
16 weight: Weighting scheme.
17 'none' (mean no weight)
18 'var' (1/var(spec) weighted)
19 'tsys' (1/Tsys**2 weighted)
20 'tint' (integration time weighted)
21 'tintsys' (Tint/Tsys**2)
22 'median' ( median averaging)
23 align: align the spectra in velocity before averaging. It takes
24 the time of the first spectrum in the first scantable
25 as reference time.
26 Example:
27 # return a time averaged scan from scana and scanb
28 # without using a mask
29 scanav = average_time(scana,scanb)
30 # or equivalent
31 # scanav = average_time([scana, scanb])
32 # return the (time) averaged scan, i.e. the average of
33 # all correlator cycles
34 scanav = average_time(scan, scanav=True)
35 """
36 scanav = False
37 if kwargs.has_key('scanav'):
38 scanav = kwargs.get('scanav')
39 weight = 'tint'
40 if kwargs.has_key('weight'):
41 weight = kwargs.get('weight')
42 mask = ()
43 if kwargs.has_key('mask'):
44 mask = kwargs.get('mask')
45 align = False
46 if kwargs.has_key('align'):
47 align = kwargs.get('align')
48 compel = False
49 if kwargs.has_key('compel'):
50 compel = kwargs.get('compel')
51 varlist = vars()
52 if isinstance(args[0],list):
53 lst = args[0]
54 elif isinstance(args[0],tuple):
55 lst = list(args[0])
56 else:
57 lst = list(args)
58
59 del varlist["kwargs"]
60 varlist["args"] = "%d scantables" % len(lst)
61 # need special formatting here for history...
62
63 from asap._asap import stmath
64 stm = stmath()
65 for s in lst:
66 if not isinstance(s,scantable):
67 msg = "Please give a list of scantables"
68 if rcParams['verbose']:
69 #print msg
70 asaplog.push(msg)
71 print_log('ERROR')
72 return
73 else:
74 raise TypeError(msg)
75 if scanav: scanav = "SCAN"
76 else: scanav = "NONE"
77 alignedlst = []
78 if align:
79 refepoch = lst[0].get_time(0)
80 for scan in lst:
81 alignedlst.append(scan.freq_align(refepoch,insitu=False))
82 else:
83 alignedlst = lst
84 if weight.upper() == 'MEDIAN':
85 # median doesn't support list of scantables - merge first
86 merged = None
87 if len(alignedlst) > 1:
88 merged = merge(alignedlst)
89 else:
90 merged = alignedlst[0]
91 s = scantable(stm._averagechannel(merged, 'MEDIAN', scanav))
92 del merged
93 else:
94 #s = scantable(stm._average(alignedlst, mask, weight.upper(), scanav))
95 s = scantable(stm._new_average(alignedlst, compel, mask, weight.upper(), scanav))
96 s._add_history("average_time",varlist)
97 print_log()
98 return s
99
100def quotient(source, reference, preserve=True):
101 """
102 Return the quotient of a 'source' (signal) scan and a 'reference' scan.
103 The reference can have just one scan, even if the signal has many. Otherwise
104 they must have the same number of scans.
105 The cursor of the output scan is set to 0
106 Parameters:
107 source: the 'on' scan
108 reference: the 'off' scan
109 preserve: you can preserve (default) the continuum or
110 remove it. The equations used are
111 preserve: Output = Toff * (on/off) - Toff
112 remove: Output = Toff * (on/off) - Ton
113 """
114 varlist = vars()
115 from asap._asap import stmath
116 stm = stmath()
117 stm._setinsitu(False)
118 s = scantable(stm._quotient(source, reference, preserve))
119 s._add_history("quotient",varlist)
120 print_log()
121 return s
122
123def dototalpower(calon, caloff, tcalval=0.0):
124 """
125 Do calibration for CAL on,off signals.
126 Adopted from GBTIDL dototalpower
127 Parameters:
128 calon: the 'cal on' subintegration
129 caloff: the 'cal off' subintegration
130 tcalval: user supplied Tcal value
131 """
132 varlist = vars()
133 from asap._asap import stmath
134 stm = stmath()
135 stm._setinsitu(False)
136 s = scantable(stm._dototalpower(calon, caloff, tcalval))
137 s._add_history("dototalpower",varlist)
138 print_log()
139 return s
140
141def dosigref(sig, ref, smooth, tsysval=0.0, tauval=0.0):
142 """
143 Calculate a quotient (sig-ref/ref * Tsys)
144 Adopted from GBTIDL dosigref
145 Parameters:
146 sig: on source data
147 ref: reference data
148 smooth: width of box car smoothing for reference
149 tsysval: user specified Tsys (scalar only)
150 tauval: user specified Tau (required if tsysval is set)
151 """
152 varlist = vars()
153 from asap._asap import stmath
154 stm = stmath()
155 stm._setinsitu(False)
156 s = scantable(stm._dosigref(sig, ref, smooth, tsysval, tauval))
157 s._add_history("dosigref",varlist)
158 print_log()
159 return s
160
161def calps(scantab, scannos, smooth=1, tsysval=0.0, tauval=0.0, tcalval=0.0):
162 """
163 Calibrate GBT position switched data
164 Adopted from GBTIDL getps
165 Currently calps identify the scans as position switched data if they
166 contain '_ps' in the source name. The data must contains 'CAL' signal
167 on/off in each integration. To identify 'CAL' on state, the word, 'calon'
168 need to be present in the source name field.
169 (GBT MS data reading process to scantable automatically append these
170 id names to the source names)
171
172 Parameters:
173 scantab: scantable
174 scannos: list of scan numbers
175 smooth: optional box smoothing order for the reference
176 (default is 1 = no smoothing)
177 tsysval: optional user specified Tsys (default is 0.0,
178 use Tsys in the data)
179 tauval: optional user specified Tau
180 tcalval: optional user specified Tcal (default is 0.0,
181 use Tcal value in the data)
182 """
183 varlist = vars()
184 # check for the appropriate data
185 s = scantab.get_scan('*_ps*')
186 if s is None:
187 msg = "The input data appear to contain no position-switch mode data."
188 if rcParams['verbose']:
189 #print msg
190 asaplog.push(msg)
191 print_log('ERROR')
192 return
193 else:
194 raise TypeError(msg)
195 ssub = s.get_scan(scannos)
196 if ssub is None:
197 msg = "No data was found with given scan numbers!"
198 if rcParams['verbose']:
199 #print msg
200 asaplog.push(msg)
201 print_log('ERROR')
202 return
203 else:
204 raise TypeError(msg)
205 ssubon = ssub.get_scan('*calon')
206 ssuboff = ssub.get_scan('*[^calon]')
207 if ssubon.nrow() != ssuboff.nrow():
208 msg = "mismatch in numbers of CAL on/off scans. Cannot calibrate. Check the scan numbers."
209 if rcParams['verbose']:
210 #print msg
211 asaplog.push(msg)
212 print_log('ERROR')
213 return
214 else:
215 raise TypeError(msg)
216 cals = dototalpower(ssubon, ssuboff, tcalval)
217 sig = cals.get_scan('*ps')
218 ref = cals.get_scan('*psr')
219 if sig.nscan() != ref.nscan():
220 msg = "mismatch in numbers of on/off scans. Cannot calibrate. Check the scan numbers."
221 if rcParams['verbose']:
222 #print msg
223 asaplog.push(msg)
224 print_log('ERROR')
225 return
226 else:
227 raise TypeError(msg)
228
229 #for user supplied Tsys
230 if tsysval>0.0:
231 if tauval<=0.0:
232 msg = "Need to supply a valid tau to use the supplied Tsys"
233 if rcParams['verbose']:
234 #print msg
235 asaplog.push(msg)
236 print_log('ERROR')
237 return
238 else:
239 raise TypeError(msg)
240 else:
241 sig.recalc_azel()
242 ref.recalc_azel()
243 #msg = "Use of user specified Tsys is not fully implemented yet."
244 #if rcParams['verbose']:
245 # print msg
246 # return
247 #else:
248 # raise TypeError(msg)
249 # use get_elevation to get elevation and
250 # calculate a scaling factor using the formula
251 # -> tsys use to dosigref
252
253 #ress = dosigref(sig, ref, smooth, tsysval)
254 ress = dosigref(sig, ref, smooth, tsysval, tauval)
255 ress._add_history("calps", varlist)
256 print_log()
257 return ress
258
259def calnod(scantab, scannos=[], smooth=1, tsysval=0.0, tauval=0.0, tcalval=0.0):
260 """
261 Do full (but a pair of scans at time) processing of GBT Nod data
262 calibration.
263 Adopted from GBTIDL's getnod
264 Parameters:
265 scantab: scantable
266 scannos: a pair of scan numbers, or the first scan number of the pair
267 smooth: box car smoothing order
268 tsysval: optional user specified Tsys value
269 tauval: optional user specified tau value (not implemented yet)
270 tcalval: optional user specified Tcal value
271 """
272 varlist = vars()
273 from asap._asap import stmath
274 stm = stmath()
275 stm._setinsitu(False)
276
277 # check for the appropriate data
278 s = scantab.get_scan('*_nod*')
279 if s is None:
280 msg = "The input data appear to contain no Nod observing mode data."
281 if rcParams['verbose']:
282 #print msg
283 asaplog.push(msg)
284 print_log('ERROR')
285 return
286 else:
287 raise TypeError(msg)
288
289 # need check correspondance of each beam with sig-ref ...
290 # check for timestamps, scan numbers, subscan id (not available in
291 # ASAP data format...). Assume 1st scan of the pair (beam 0 - sig
292 # and beam 1 - ref...)
293 # First scan number of paired scans or list of pairs of
294 # scan numbers (has to have even number of pairs.)
295
296 #data splitting
297 scan1no = scan2no = 0
298
299 if len(scannos)==1:
300 scan1no = scannos[0]
301 scan2no = scannos[0]+1
302 pairScans = [scan1no, scan2no]
303 else:
304 #if len(scannos)>2:
305 # msg = "calnod can only process a pair of nod scans at time."
306 # if rcParams['verbose']:
307 # print msg
308 # return
309 # else:
310 # raise TypeError(msg)
311 #
312 #if len(scannos)==2:
313 # scan1no = scannos[0]
314 # scan2no = scannos[1]
315 pairScans = list(scannos)
316
317 if tsysval>0.0:
318 if tauval<=0.0:
319 msg = "Need to supply a valid tau to use the supplied Tsys"
320 if rcParams['verbose']:
321 #print msg
322 asaplog.push(msg)
323 print_log('ERROR')
324 return
325 else:
326 raise TypeError(msg)
327 else:
328 scantab.recalc_azel()
329 resspec = scantable(stm._donod(scantab, pairScans, smooth, tsysval,tauval,tcalval))
330 resspec._add_history("calnod",varlist)
331 print_log()
332 return resspec
333
334def calfs(scantab, scannos=[], smooth=1, tsysval=0.0, tauval=0.0, tcalval=0.0):
335 """
336 Calibrate GBT frequency switched data.
337 Adopted from GBTIDL getfs.
338 Currently calfs identify the scans as frequency switched data if they
339 contain '_fs' in the source name. The data must contains 'CAL' signal
340 on/off in each integration. To identify 'CAL' on state, the word, 'calon'
341 need to be present in the source name field.
342 (GBT MS data reading via scantable automatically append these
343 id names to the source names)
344
345 Parameters:
346 scantab: scantable
347 scannos: list of scan numbers
348 smooth: optional box smoothing order for the reference
349 (default is 1 = no smoothing)
350 tsysval: optional user specified Tsys (default is 0.0,
351 use Tsys in the data)
352 tauval: optional user specified Tau
353 """
354 varlist = vars()
355 from asap._asap import stmath
356 stm = stmath()
357 stm._setinsitu(False)
358
359# check = scantab.get_scan('*_fs*')
360# if check is None:
361# msg = "The input data appear to contain no Nod observing mode data."
362# if rcParams['verbose']:
363# print msg
364# return
365# else:
366# raise TypeError(msg)
367 s = scantab.get_scan(scannos)
368 del scantab
369
370 resspec = scantable(stm._dofs(s, scannos, smooth, tsysval,tauval,tcalval))
371 resspec._add_history("calfs",varlist)
372 print_log()
373 return resspec
374
375def simple_math(left, right, op='add', tsys=True):
376 """
377 Apply simple mathematical binary operations to two
378 scan tables, returning the result in a new scan table.
379 The operation is applied to both the correlations and the TSys data
380 The cursor of the output scan is set to 0
381 Parameters:
382 left: the 'left' scan
383 right: the 'right' scan
384 op: the operation: 'add' (default), 'sub', 'mul', 'div'
385 tsys: if True (default) then apply the operation to Tsys
386 as well as the data
387 """
388 #print "simple_math is deprecated use +=/* instead."
389 asaplog.push( "simple_math is deprecated use +=/* instead." )
390 print_log('WARN')
391
392def merge(*args):
393 """
394 Merge a list of scanatables, or comma-sperated scantables into one
395 scnatble.
396 Parameters:
397 A list [scan1, scan2] or scan1, scan2.
398 Example:
399 myscans = [scan1, scan2]
400 allscans = merge(myscans)
401 # or equivalent
402 sameallscans = merge(scan1, scan2)
403 """
404 varlist = vars()
405 if isinstance(args[0],list):
406 lst = tuple(args[0])
407 elif isinstance(args[0],tuple):
408 lst = args[0]
409 else:
410 lst = tuple(args)
411 varlist["args"] = "%d scantables" % len(lst)
412 # need special formatting her for history...
413 from asap._asap import stmath
414 stm = stmath()
415 for s in lst:
416 if not isinstance(s,scantable):
417 msg = "Please give a list of scantables"
418 if rcParams['verbose']:
419 #print msg
420 asaplog.push(msg)
421 print_log('ERROR')
422 return
423 else:
424 raise TypeError(msg)
425 s = scantable(stm._merge(lst))
426 s._add_history("merge", varlist)
427 print_log()
428 return s
429
Note: See TracBrowser for help on using the repository browser.