casa: predictcomp_helper.py Source File

Go to the documentation of this file.
00001 import pylab as pl
00002 from taskinit import *
00003 import shutil
00004 from setjy_helper import *
00005 
00006 def predictSolarObjectCompList(objname, epoch, freqs, prefix):
00007     """
00008     predictcomp functionality using the new models
00009     set flux density of a solar system object using Bryan Butler's new
00010     python model calculation code.
00011     """
00012     retval = True
00013     cleanupcomps = False # leave genenerated cl files
00014     nfreqs=-1
00015 
00016     (myms, mytb, mycl,myme) = gentools(['ms','tb','cl','me'])
00017     #freqinc=freqs[0]*1e-6
00018     if len(freqs) == 1:
00019       freqinc=freqs[0]*1e-4
00020     else:
00021       freqinc=abs(freqs[1]-freqs[0])*1e-2
00022     freqlist=[]
00023     for freq in freqs:
00024       minf = freq-freqinc
00025       maxf = freq+freqinc
00026       freqlist.append([minf,maxf]) 
00027     #mepoch = myme.epoch('UTC', epoch)
00028     if epoch['m0']['value']==0.0:
00029         casalog.post('Invalid epoch, '+str(epoch['m0']['value'])+str(epoch['m0']['unit']),'SEVERE');
00030         raise Exception, "Error"
00031     epochv = epoch['m0']['value'] 
00032 
00033 
00034     # turn user input epoch to mjd
00035 
00036     import solar_system_setjy as ss_setjy
00037     #print "sending objname=",objname, " epochv=",epochv, " freqlist=",freqlist
00038     observatory='ALMA'
00039     (errcodes, fluxes, fluxerrs, sizes, dirs)=\
00040        ss_setjy.solar_system_fd(source_name=objname, MJDs=[epochv], frequencies=freqlist, observatory=observatory, casalog=casalog)
00041   
00042     #print "fluxes from ss_setjy=",fluxes
00043     #if errcodes[0][0] > 0:
00044     #    raise ValueError("cannot determined flux")
00045     
00046     reterr=testerrs(errcodes[0],objname) 
00047     if reterr == 2: 
00048         #raise Exception, "Error" 
00049         casalog.post("Flux densities cannot be determined","SEVERE")
00050         raise Exception
00051 
00052     dirstring = [dirs[0]['refer'],qa.tos(dirs[0]['m0']),qa.tos(dirs[0]['m1'])]
00053     # setup componentlists
00054     # need to set per dir
00055     # if scalebychan=F, len(freqs) corresponds to nspw selected
00056 
00057     #clpath='/tmp/'
00058     clpath='./'
00059     # use the first input frequency
00060     freqlabel = '%.3fGHz' % (freqs[0]/1.e9)
00061     tmlabel = '%.1fd' % epoch['m0']['value']
00062     clabel = objname+'_spw0_'+freqlabel+'_'+tmlabel
00063     clname = clabel+'.cl'
00064 
00065     if(os.path.exists(clname)):
00066         print "Removing previous cl file,", clname
00067         try:
00068             shutil.rmtree(clname)
00069         except:
00070             print "shutil.rmtree failed" 
00071     index= 2.0
00072     sptype = 'spectral index'
00073     #index= 0.0
00074     #sptype = 'constant'
00075     
00076     #print "fluxes=",fluxes
00077     #print "fluxerrs-=",fluxerrs
00078     #print "sizes=",sizes
00079     #print "dirs=",dirs 
00080     mycl.addcomponent(flux=fluxes[0][0],fluxunit='Jy', polarization="Stokes", dir=dirs[0],
00081          shape='disk', majoraxis=str(sizes[0][0])+'arcsec', minoraxis=str(sizes[0][1])+'arcsec',
00082          positionangle=str(sizes[0][2])+'arcsec', freq=['LSRK',str(freqs[0])+'Hz'],
00083          spectrumtype=sptype, index=index, label=clabel)
00084     # set tabular, use original input freqs
00085     if len(freqs)>1:
00086         mycl.setspectrum(which=0, type='tabular', tabularfreqs=freqs, tabularflux=fluxes[0],
00087                            tabularframe='TOPO')
00088 
00089     mycl.rename(clname)
00090     mycl.close(False)
00091     return (clname)