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.


File:
1 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')
Note: See TracChangeset for help on using the changeset viewer.