source: branches/asap-3.x/python/asapmath.py@ 3150

Last change on this file since 3150 was 1591, checked in by Malte Marquarding, 15 years ago

Removed deprecated function

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