source: branches/mergetest/python/asapmath.py @ 1777

Last change on this file since 1777 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.