Changeset 1631 for branches/alma/python


Ignore:
Timestamp:
09/08/09 18:47:19 (15 years ago)
Author:
Takeshi Nakazato
Message:

New Development: No

JIRA Issue: Yes CAS-1424

Ready to Release: Yes

Interface Changes: Yes

What Interface Changed: Added 'verify' parameter for asapmath.calps, asapmath.calnod, asapmath.calfs, and scantable.auto_quotient

Test Programs: do calibration using sdaverage with verify=True

Put in Release Notes: Yes

Module(s): task sdaverage and sdcal

Description: Describe your changes here...

I have added verify parameter to verify calibration result.
If verify=True, methods plot pre- and post-calibration spectra and asks user
to be acceptable calibration or not.
If not, calibration will not be applied.


Location:
branches/alma/python
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • branches/alma/python/asapmath.py

    r1614 r1631  
    44from asap import selector
    55from asap import asaplog
     6from asap import asaplotgui
    67
    78def average_time(*args, **kwargs):
    … …  
    159160    return s
    160161
    161 def calps(scantab, scannos, smooth=1, tsysval=0.0, tauval=0.0, tcalval=0.0):
     162def calps(scantab, scannos, smooth=1, tsysval=0.0, tauval=0.0, tcalval=0.0, verify=False):
    162163    """
    163164    Calibrate GBT position switched data
    … …  
    253254    #ress = dosigref(sig, ref, smooth, tsysval)
    254255    ress = dosigref(sig, ref, smooth, tsysval, tauval)
     256    ###
     257    if verify:
     258        # get data
     259        import numpy
     260        precal={}
     261        postcal=[]
     262        keys=['ps','ps_calon','psr','psr_calon']
     263        ifnos=list(ssub.getifnos())
     264        polnos=list(ssub.getpolnos())
     265        sel=selector()
     266        for i in range(2):
     267            ss=ssuboff.get_scan('*'+keys[2*i])
     268            ll=[]
     269            for j in range(len(ifnos)):
     270                for k in range(len(polnos)):
     271                    sel.set_ifs(ifnos[j])
     272                    sel.set_polarizations(polnos[k])
     273                    try:
     274                        ss.set_selection(sel)
     275                    except:
     276                        continue
     277                    ll.append(numpy.array(ss._getspectrum(0)))
     278                    sel.reset()
     279                    ss.set_selection()
     280            precal[keys[2*i]]=ll
     281            del ss
     282            ss=ssubon.get_scan('*'+keys[2*i+1])
     283            ll=[]
     284            for j in range(len(ifnos)):
     285                for k in range(len(polnos)):
     286                    sel.set_ifs(ifnos[j])
     287                    sel.set_polarizations(polnos[k])
     288                    try:
     289                        ss.set_selection(sel)
     290                    except:
     291                        continue
     292                    ll.append(numpy.array(ss._getspectrum(0)))
     293                    sel.reset()
     294                    ss.set_selection()
     295            precal[keys[2*i+1]]=ll
     296            del ss
     297        for j in range(len(ifnos)):
     298            for k in range(len(polnos)):
     299                sel.set_ifs(ifnos[j])
     300                sel.set_polarizations(polnos[k])
     301                try:
     302                    ress.set_selection(sel)
     303                except:
     304                    continue
     305                postcal.append(numpy.array(ress._getspectrum(0)))
     306                sel.reset()
     307                ress.set_selection()
     308        del sel
     309        # plot
     310        print_log()
     311        asaplog.push('Plot only first spectrum for each [if,pol] pairs to verify calibration.')
     312        print_log('WARN')
     313        p=asaplotgui.asaplotgui()
     314        #nr=min(6,len(ifnos)*len(polnos))
     315        nr=len(ifnos)*len(polnos)
     316        titles=[]
     317        btics=[]
     318        if nr<4:
     319            p.set_panels(rows=nr,cols=2,nplots=2*nr,ganged=False)
     320            for i in range(2*nr):
     321                b=False
     322                if i >= 2*nr-2:
     323                    b=True
     324                btics.append(b)
     325        elif nr==4:
     326            p.set_panels(rows=2,cols=4,nplots=8,ganged=False)
     327            for i in range(2*nr):
     328                b=False
     329                if i >= 2*nr-4:
     330                    b=True
     331                btics.append(b)
     332        elif nr<7:
     333            p.set_panels(rows=3,cols=4,nplots=2*nr,ganged=False)
     334            for i in range(2*nr):
     335                if i >= 2*nr-4:
     336                    b=True
     337                btics.append(b)
     338        else:
     339            print_log()
     340            asaplog.push('Only first 6 [if,pol] pairs are plotted.')
     341            print_log('WARN')
     342            nr=6
     343            for i in range(2*nr):
     344                b=False
     345                if i >= 2*nr-4:
     346                    b=True
     347                btics.append(b)
     348            p.set_panels(rows=3,cols=4,nplots=2*nr,ganged=False)
     349        for i in range(nr):
     350            p.subplot(2*i)
     351            p.color=0
     352            title='raw data IF%s POL%s' % (ifnos[int(i/len(polnos))],polnos[i%len(polnos)])
     353            titles.append(title)
     354            #p.set_axes('title',title,fontsize=40)
     355            ymin=1.0e100
     356            ymax=-1.0e100
     357            nchan=s.nchan()
     358            edge=int(nchan*0.01)
     359            for j in range(4):
     360                spmin=min(precal[keys[j]][i][edge:nchan-edge])
     361                spmax=max(precal[keys[j]][i][edge:nchan-edge])
     362                ymin=min(ymin,spmin)
     363                ymax=max(ymax,spmax)
     364            for j in range(4):
     365                if i==0:
     366                    p.set_line(label=keys[j])
     367                else:
     368                    p.legend()
     369                p.plot(precal[keys[j]][i])
     370            p.axes.set_ylim(ymin-0.1*abs(ymin),ymax+0.1*abs(ymax))
     371            if not btics[2*i]:
     372                p.axes.set_xticks([])
     373            p.subplot(2*i+1)
     374            p.color=0
     375            title='cal data IF%s POL%s' % (ifnos[int(i/len(polnos))],polnos[i%len(polnos)])
     376            titles.append(title)
     377            #p.set_axes('title',title)
     378            p.legend()
     379            ymin=postcal[i][edge:nchan-edge].min()
     380            ymax=postcal[i][edge:nchan-edge].max()
     381            p.plot(postcal[i])
     382            p.axes.set_ylim(ymin-0.1*abs(ymin),ymax+0.1*abs(ymax))
     383            if not btics[2*i+1]:
     384                p.axes.set_xticks([])
     385        for i in range(2*nr):
     386            p.subplot(i)
     387            p.set_axes('title',titles[i],fontsize='medium')
     388        x=raw_input('Accept calibration ([y]/n): ' )
     389        if x.upper() == 'N':
     390            p.unmap()
     391            del p
     392            return scabtab
     393        p.unmap()
     394        del p
     395    ###
    255396    ress._add_history("calps", varlist)
    256397    print_log()
    257398    return ress
    258399
    259 def calnod(scantab, scannos=[], smooth=1, tsysval=0.0, tauval=0.0, tcalval=0.0):
     400def calnod(scantab, scannos=[], smooth=1, tsysval=0.0, tauval=0.0, tcalval=0.0, verify=False):
    260401    """
    261402    Do full (but a pair of scans at time) processing of GBT Nod data
    … …  
    328469            scantab.recalc_azel()
    329470    resspec = scantable(stm._donod(scantab, pairScans, smooth, tsysval,tauval,tcalval))
     471    ###
     472    if verify:
     473        # get data
     474        import numpy
     475        precal={}
     476        postcal=[]
     477        keys=['nod','nod_calon']
     478        ifnos=list(scantab.getifnos())
     479        polnos=list(scantab.getpolnos())
     480        sel=selector()
     481        for i in range(2):
     482            ss=scantab.get_scan('*'+keys[i])
     483            ll=[]
     484            ll2=[]
     485            for j in range(len(ifnos)):
     486                for k in range(len(polnos)):
     487                    sel.set_ifs(ifnos[j])
     488                    sel.set_polarizations(polnos[k])
     489                    sel.set_scans(pairScans[0])
     490                    try:
     491                        ss.set_selection(sel)
     492                    except:
     493                        continue
     494                    ll.append(numpy.array(ss._getspectrum(0)))
     495                    sel.reset()
     496                    ss.set_selection()
     497                    sel.set_ifs(ifnos[j])
     498                    sel.set_polarizations(polnos[k])
     499                    sel.set_scans(pairScans[1])
     500                    try:
     501                        ss.set_selection(sel)
     502                    except:
     503                        ll.pop()
     504                        continue
     505                    ll2.append(numpy.array(ss._getspectrum(0)))
     506                    sel.reset()
     507                    ss.set_selection()
     508            key='%s%s' %(pairScans[0],keys[i].lstrip('nod'))
     509            precal[key]=ll
     510            key='%s%s' %(pairScans[1],keys[i].lstrip('nod'))
     511            precal[key]=ll2
     512            del ss
     513        keys=precal.keys()
     514        for j in range(len(ifnos)):
     515            for k in range(len(polnos)):
     516                sel.set_ifs(ifnos[j])
     517                sel.set_polarizations(polnos[k])
     518                sel.set_scans(pairScans[0])
     519                try:
     520                    resspec.set_selection(sel)
     521                except:
     522                    continue
     523                postcal.append(numpy.array(resspec._getspectrum(0)))
     524                sel.reset()
     525                resspec.set_selection()
     526        del sel
     527        # plot
     528        print_log()
     529        asaplog.push('Plot only first spectrum for each [if,pol] pairs to verify calibration.')
     530        print_log('WARN')
     531        p=asaplotgui.asaplotgui()
     532        #nr=min(6,len(ifnos)*len(polnos))
     533        nr=len(ifnos)*len(polnos)
     534        titles=[]
     535        btics=[]
     536        if nr<4:
     537            p.set_panels(rows=nr,cols=2,nplots=2*nr,ganged=False)
     538            for i in range(2*nr):
     539                b=False
     540                if i >= 2*nr-2:
     541                    b=True
     542                btics.append(b)
     543        elif nr==4:
     544            p.set_panels(rows=2,cols=4,nplots=8,ganged=False)
     545            for i in range(2*nr):
     546                b=False
     547                if i >= 2*nr-4:
     548                    b=True
     549                btics.append(b)
     550        elif nr<7:
     551            p.set_panels(rows=3,cols=4,nplots=2*nr,ganged=False)
     552            for i in range(2*nr):
     553                if i >= 2*nr-4:
     554                    b=True
     555                btics.append(b)
     556        else:
     557            print_log()
     558            asaplog.push('Only first 6 [if,pol] pairs are plotted.')
     559            print_log('WARN')
     560            nr=6
     561            for i in range(2*nr):
     562                b=False
     563                if i >= 2*nr-4:
     564                    b=True
     565                btics.append(b)
     566            p.set_panels(rows=3,cols=4,nplots=2*nr,ganged=False)
     567        for i in range(nr):
     568            p.subplot(2*i)
     569            p.color=0
     570            title='raw data IF%s POL%s' % (ifnos[int(i/len(polnos))],polnos[i%len(polnos)])
     571            titles.append(title)
     572            #p.set_axes('title',title,fontsize=40)
     573            ymin=1.0e100
     574            ymax=-1.0e100
     575            nchan=scantab.nchan()
     576            edge=int(nchan*0.01)
     577            for j in range(4):
     578                spmin=min(precal[keys[j]][i][edge:nchan-edge])
     579                spmax=max(precal[keys[j]][i][edge:nchan-edge])
     580                ymin=min(ymin,spmin)
     581                ymax=max(ymax,spmax)
     582            for j in range(4):
     583                if i==0:
     584                    p.set_line(label=keys[j])
     585                else:
     586                    p.legend()
     587                p.plot(precal[keys[j]][i])
     588            p.axes.set_ylim(ymin-0.1*abs(ymin),ymax+0.1*abs(ymax))
     589            if not btics[2*i]:
     590                p.axes.set_xticks([])
     591            p.subplot(2*i+1)
     592            p.color=0
     593            title='cal data IF%s POL%s' % (ifnos[int(i/len(polnos))],polnos[i%len(polnos)])
     594            titles.append(title)
     595            #p.set_axes('title',title)
     596            p.legend()
     597            ymin=postcal[i][edge:nchan-edge].min()
     598            ymax=postcal[i][edge:nchan-edge].max()
     599            p.plot(postcal[i])
     600            p.axes.set_ylim(ymin-0.1*abs(ymin),ymax+0.1*abs(ymax))
     601            if not btics[2*i+1]:
     602                p.axes.set_xticks([])
     603        for i in range(2*nr):
     604            p.subplot(i)
     605            p.set_axes('title',titles[i],fontsize='medium')
     606        x=raw_input('Accept calibration ([y]/n): ' )
     607        if x.upper() == 'N':
     608            p.unmap()
     609            del p
     610            return scabtab
     611        p.unmap()
     612        del p
     613    ###
    330614    resspec._add_history("calnod",varlist)
    331615    print_log()
    332616    return resspec
    333617
    334 def calfs(scantab, scannos=[], smooth=1, tsysval=0.0, tauval=0.0, tcalval=0.0):
     618def calfs(scantab, scannos=[], smooth=1, tsysval=0.0, tauval=0.0, tcalval=0.0, verify=False):
    335619    """
    336620    Calibrate GBT frequency switched data.
    … …  
    369653
    370654    resspec = scantable(stm._dofs(s, scannos, smooth, tsysval,tauval,tcalval))
     655    ###
     656    if verify:
     657        # get data
     658        ssub = s.get_scan(scannos)
     659        ssubon = ssub.get_scan('*calon')
     660        ssuboff = ssub.get_scan('*[^calon]')
     661        import numpy
     662        precal={}
     663        postcal=[]
     664        keys=['fs','fs_calon','fsr','fsr_calon']
     665        ifnos=list(ssub.getifnos())
     666        polnos=list(ssub.getpolnos())
     667        sel=selector()
     668        for i in range(2):
     669            ss=ssuboff.get_scan('*'+keys[2*i])
     670            ll=[]
     671            for j in range(len(ifnos)):
     672                for k in range(len(polnos)):
     673                    sel.set_ifs(ifnos[j])
     674                    sel.set_polarizations(polnos[k])
     675                    try:
     676                        ss.set_selection(sel)
     677                    except:
     678                        continue
     679                    ll.append(numpy.array(ss._getspectrum(0)))
     680                    sel.reset()
     681                    ss.set_selection()
     682            precal[keys[2*i]]=ll
     683            del ss
     684            ss=ssubon.get_scan('*'+keys[2*i+1])
     685            ll=[]
     686            for j in range(len(ifnos)):
     687                for k in range(len(polnos)):
     688                    sel.set_ifs(ifnos[j])
     689                    sel.set_polarizations(polnos[k])
     690                    try:
     691                        ss.set_selection(sel)
     692                    except:
     693                        continue
     694                    ll.append(numpy.array(ss._getspectrum(0)))
     695                    sel.reset()
     696                    ss.set_selection()
     697            precal[keys[2*i+1]]=ll
     698            del ss
     699        sig=resspec.get_scan('*_fs')
     700        ref=resspec.get_scan('*_fsr')
     701        for k in range(len(polnos)):
     702            for j in range(len(ifnos)):
     703                sel.set_ifs(ifnos[j])
     704                sel.set_polarizations(polnos[k])
     705                try:
     706                    sig.set_selection(sel)
     707                    postcal.append(numpy.array(sig._getspectrum(0)))
     708                except:
     709                    ref.set_selection(sel)
     710                    postcal.append(numpy.array(ref._getspectrum(0)))
     711                sel.reset()
     712                resspec.set_selection()
     713        del sel
     714        # plot
     715        print_log()
     716        asaplog.push('Plot only first spectrum for each [if,pol] pairs to verify calibration.')
     717        print_log('WARN')
     718        p=asaplotgui.asaplotgui()
     719        #nr=min(6,len(ifnos)*len(polnos))
     720        nr=len(ifnos)/2*len(polnos)
     721        titles=[]
     722        btics=[]
     723        if nr>3:
     724            print_log()
     725            asaplog.push('Only first 3 [if,pol] pairs are plotted.')
     726            print_log('WARN')
     727            nr=3
     728        p.set_panels(rows=nr,cols=3,nplots=3*nr,ganged=False)
     729        for i in range(3*nr):
     730            b=False
     731            if i >= 3*nr-3:
     732                b=True
     733            btics.append(b)
     734        for i in range(nr):
     735            p.subplot(3*i)
     736            p.color=0
     737            title='raw data IF%s,%s POL%s' % (ifnos[2*int(i/len(polnos))],ifnos[2*int(i/len(polnos))+1],polnos[i%len(polnos)])
     738            titles.append(title)
     739            #p.set_axes('title',title,fontsize=40)
     740            ymin=1.0e100
     741            ymax=-1.0e100
     742            nchan=s.nchan()
     743            edge=int(nchan*0.01)
     744            for j in range(4):
     745                spmin=min(precal[keys[j]][i][edge:nchan-edge])
     746                spmax=max(precal[keys[j]][i][edge:nchan-edge])
     747                ymin=min(ymin,spmin)
     748                ymax=max(ymax,spmax)
     749            for j in range(4):
     750                if i==0:
     751                    p.set_line(label=keys[j])
     752                else:
     753                    p.legend()
     754                p.plot(precal[keys[j]][i])
     755            p.axes.set_ylim(ymin-0.1*abs(ymin),ymax+0.1*abs(ymax))
     756            if not btics[3*i]:
     757                p.axes.set_xticks([])
     758            p.subplot(3*i+1)
     759            p.color=0
     760            title='sig data IF%s POL%s' % (ifnos[2*int(i/len(polnos))],polnos[i%len(polnos)])
     761            titles.append(title)
     762            #p.set_axes('title',title)
     763            p.legend()
     764            ymin=postcal[2*i][edge:nchan-edge].min()
     765            ymax=postcal[2*i][edge:nchan-edge].max()
     766            p.plot(postcal[2*i])
     767            p.axes.set_ylim(ymin-0.1*abs(ymin),ymax+0.1*abs(ymax))
     768            if not btics[3*i+1]:
     769                p.axes.set_xticks([])
     770            p.subplot(3*i+2)
     771            p.color=0
     772            title='ref data IF%s POL%s' % (ifnos[2*int(i/len(polnos))+1],polnos[i%len(polnos)])
     773            titles.append(title)
     774            #p.set_axes('title',title)
     775            p.legend()
     776            ymin=postcal[2*i+1][edge:nchan-edge].min()
     777            ymax=postcal[2*i+1][edge:nchan-edge].max()
     778            p.plot(postcal[2*i+1])
     779            p.axes.set_ylim(ymin-0.1*abs(ymin),ymax+0.1*abs(ymax))
     780            if not btics[3*i+2]:
     781                p.axes.set_xticks([])
     782        for i in range(3*nr):
     783            p.subplot(i)
     784            p.set_axes('title',titles[i],fontsize='medium')
     785        x=raw_input('Accept calibration ([y]/n): ' )
     786        if x.upper() == 'N':
     787            p.unmap()
     788            del p
     789            return scabtab
     790        p.unmap()
     791        del p
     792    ###
    371793    resspec._add_history("calfs",varlist)
    372794    print_log()
    … …  
    428850    return s
    429851
     852## def apexcal( scantab, calmode ):
     853##     """
     854##     Calibrate APEX data.
     855   
     856##     Parameters:
     857##         scantab:       scantable
     858##         calmode:       calibration mode
     859##     """
     860##     from asap._asap import stmath
     861##     stm = stmath()
     862##     if ( calmode == 'ps' ):
     863##         asaplog.push( 'APEX position-switch calibration' )
     864##         print_log()
     865##     elif ( calmode == 'fs' ):
     866##         asaplog.push( 'APEX frequency-switch calibration' )
     867##         print_log()
     868##     elif ( calmode == 'wob' ):
     869##         asaplog.push( 'APEX wobbler-switch calibration' )
     870##         print_log()
     871##     elif ( calmode == 'otf' ):
     872##         asaplog.push( 'APEX On-The-Fly calibration' )
     873##         print_log()
     874##     else:
     875##         asaplog.push( '%s: unknown calibration mode' % calmode )
     876##         print_log('ERROR')

  • branches/alma/python/scantable.py

    r1629 r1631  
    18161816        s._add_history("set_sourcetype", varlist)
    18171817
    1818     def auto_quotient(self, preserve=True, mode='paired'):
     1818    def auto_quotient(self, preserve=True, mode='paired', verify=False):
    18191819        """
    18201820        This function allows to build quotients automatically.
    1821         It assumes the observation to have the same numer of
     1821        It assumes the observation to have the same number of
    18221822        "ons" and "offs"
    18231823        Parameters:
Note: See TracChangeset for help on using the changeset viewer.