casa: task_mosaic.py Source File

Go to the documentation of this file.
00001 import os
00002 from taskinit import *
00003 
00004 def mosaic(vis,imagename,mode,alg,imsize,cell,phasecenter,stokes,niter,gain,threshold,mask,cleanbox,nchan,start,width,field,spw,timerange,restfreq,sdimage,modelimage,weighting,mosweight,rmode,robust,ftmachine,cyclefactor,cyclespeedup,scaletype,minpb,sigma,targetflux,constrainflux,prior,negcomponent,scales,npercycle,npixels, noise):
00005         """Calculate a multi-field deconvolved image with selected algorithm:
00006 
00007         """
00008 
00009         ###
00010         #Python script
00011         casalog.origin('mosaic')
00012 
00013         if(mask==[]):
00014                 mask=['']
00015         isInteractive=False;
00016         alg = 'mf'+ alg
00017         if (cleanbox== 'interactive'):
00018                 isInteractive=True
00019                 if( (mask=='') or (mask==[''])):
00020                         mask=[imagename+'_interactive.mask']
00021 
00022         #parameter_printvalues(arg_names,arg_values,arg_types)
00023 
00024         try:
00025 ###Set image parameters
00026                 imMos=imtool()
00027                 casalog.origin('mosaic')
00028                 if((len(imagename)==0) or (imagename.isspace())):
00029                         raise Exception, 'Cannot proceed with blank imagename'
00030                 if (field==''):
00031                         field='*'
00032 
00033                 if(type(phasecenter)==str):
00034                         ### blank means take field 0
00035                         if (phasecenter==''):
00036                                 phasecenter=ms.msseltoindex(vis,field=field)['field'][0]
00037                         else:
00038                                 tmppc=phasecenter
00039                                 try:
00040                                         if(len(ms.msseltoindex(vis, field=phasecenter)['field']) > 0):
00041                                                 tmppc=ms.msseltoindex(vis, field=phasecenter)['field'][0]
00042                                         ##succesful must be string like '0' or 'NGC*'
00043                                 except Exception, instance:
00044                                         ##failed must be a string 'J2000 18h00m00 10d00m00'
00045                                         tmppc=phasecenter
00046                                 phasecenter=tmppc
00047 
00048                 spwindex=-1;
00049                 if( (spw==-1) or (spw=='-1')  or (spw=='') or (spw==' ')):
00050                         spwindex=-1
00051                 else:
00052                         spwindex=ms.msseltoindex(vis, spw=spw)['spw'].tolist()
00053 
00054                 fieldindex=ms.msseltoindex(vis,field=field)['field'].tolist()
00055                 if ((type(vis)==str) & (os.path.exists(vis))):
00056                         imMos.open(vis, usescratch=mosweight)
00057                 else:
00058                         raise Exception, 'Visibility data set not found - please verify the name'
00059                 imMos.setvp(dovp=True)
00060                 if ((type(imsize)==int)):
00061                         imsize=[imsize,imsize]
00062                 if ((type(cell)==int) | (type(cell)==float) | (type(cell)==str)):
00063                         cell=[cell,cell]
00064                 if ((len(imsize)==1)): 
00065                         imsize=[imsize[0],imsize[0]]
00066                 if ((len(cell)==1)):
00067                         cell=[cell[0],cell[0]]
00068                 
00069                 cellx=cell[0]
00070                 celly=cell[1]
00071                 if((type(cell[0])==int) or (type(cell[0])==float)):
00072                         cellx=qa.quantity(cell[0], 'arcsec')
00073                         celly=qa.quantity(cell[1], 'arcsec')
00074 
00075 
00076                 print 'phasecenter= ', phasecenter
00077                 #parameter_printvalues(arg_names,arg_values,arg_types)
00078                 if(mode=='frequency'):
00079                         ##check that start and step have units
00080                         if(qa.quantity(start)['unit'].find('Hz') < 1):
00081                                 raise TypeError, "start parameter is not a valid frequency quantity "
00082                         if(qa.quantity(width)['unit'].find('Hz') < 1):
00083                                 raise TypeError, "start parameter is not a valid frequency quantity "   
00084                 elif(mode=='velocity'):
00085                         ##check that start and step have units
00086                         if(qa.quantity(start)['unit'].find('m/s') < 0):
00087                                 raise TypeError, "start parameter is not a valid velocity quantity "
00088                         if(qa.quantity(width)['unit'].find('m/s') < 0):
00089                                 raise TypeError, "start parameter is not a valid velocity quantity "
00090                 else:
00091                         if((type(width) != int) 
00092                            or (type(start) != int)):
00093                                 raise TypeError, "start, width have to be  integers with mode %s" %mode
00094 
00095                 imMos.defineimage(nx=imsize[0],ny=imsize[1],cellx=cellx,celly=celly,mode=mode,nchan=nchan,start=start,step=width,phasecenter=phasecenter,spw=spwindex,stokes=stokes, restfreq=restfreq)
00096 
00097 
00098 ###try to minimize the number of data selection
00099                 ### when natural mosweight is irrelevant                
00100                 if(weighting=='natural'):
00101                         mosweight=False 
00102                 imMos.selectvis(field=field,spw=spw,time=timerange, usescratch=True)
00103                 imMos.weight(type=weighting,rmode=rmode,robust=robust,npixels=npixels, noise=noise, mosaic=mosweight)
00104 
00105 
00106 #####
00107 ###Set clean region (file or list)
00108                 if (((mask==['']) or (mask=='')) & ( (len(cleanbox) !=0) and (type(cleanbox)==list))):
00109                         setcleanbox=False
00110                         if (type(cleanbox[0])==int):
00111                                 cleanbox=[cleanbox]
00112                                 setcleanbox = True
00113                         if ((type(cleanbox[0])==str) & (cleanbox[0]!='') & (cleanbox != 'interactive')):
00114                                 cleanbox=readboxfile(cleanbox)
00115                                 setcleanbox = True
00116 
00117                         if setcleanbox :
00118                                 mask=[imagename+'.cleanbox.mask']
00119                                 if(os.path.exists(mask[0])):
00120                                         print 'mask',mask[0], 'exists on disk, will add new mask to it'
00121                                         casalog.post('mask '+mask[0]+' exists on disk, will add new mask to it', 'WARN')
00122                                 imMos.regionmask(mask=mask[0], boxes=cleanbox)
00123 
00124 #####
00125                 if (ftmachine!='ft'): imMos.setoptions(ftmachine=ftmachine,padding=1.0)
00126                 #imMos.setmfcontrol(scaletype=scaletype,minpb=minpb)
00127                 #hardwire minpb for now
00128                 sclt=scaletype
00129                 if(scaletype=='PBCOR'):
00130                         sclt='NONE'
00131                 imMos.setmfcontrol(stoplargenegatives=negcomponent,scaletype=sclt,minpb=minpb,cyclefactor=cyclefactor,cyclespeedup=cyclespeedup,fluxscale=[imagename+'.flux'])
00132 
00133                 if (modelimage==''): modelimage=imagename+'.model'
00134 
00135 #               imMos.make(modelimage)
00136 # Using SD image as model
00137                 if (sdimage!=''): 
00138                         if (mask[0]==''):  
00139                                 imMos.makemodelfromsd(sdimage=sdimage,modelimage=modelimage,maskimage=imagename+'.mask')
00140                                 mask[0]=imagename+'.mask'
00141                         else:
00142                                 imMos.makemodelfromsd(sdimage=sdimage,modelimage=modelimage,maskimage=imagename+'.sdmask')
00143                                 mask[0]=imagename+'.sdmask'
00144 #####
00145 ###Determine algorithm to use
00146                 if (alg=='mfentropy' or alg=='mfemptiness'):
00147                         imMos.mem(algorithm=alg,niter=niter,sigma=sigma,targetflux=targetflux,constrainflux=constrainflux,keepfixed=[False],prior=prior,model=[modelimage],residual=[imagename+'.residual'],image=[imagename+'.image'],mask=mask)
00148                 elif (alg=='mfmultiscale'):
00149                         if (len(scales)==1):
00150                                 imMos.setscales('nscales',nscales=scales)
00151                         else:
00152                                 imMos.setscales('uservector',uservector=scales)
00153 
00154                         imMos.clean(algorithm=alg,niter=niter,gain=gain,threshold=qa.quantity(threshold,'mJy'),model=[modelimage],residual=[imagename+'.residual'],image=[imagename+'.image'],mask=mask,interactive=isInteractive,npercycle=npercycle)
00155                 else:
00156                         imMos.clean(algorithm=alg,niter=niter,gain=gain,threshold=qa.quantity(threshold,'mJy'),model=[modelimage],residual=[imagename+'.residual'],image=[imagename+'.image'],mask=mask,interactive=isInteractive,npercycle=npercycle)
00157 
00158                 imMos.close()
00159                 del imMos
00160                 #if (scaletype=='PBCOR'):
00161                 #       print 'dividing by fluxscale image',imagename+'.flux'
00162                 #       ia.open(imagename+'.image')
00163                 #       result=imagename+'.image'
00164                 #       fluxscale_image=imagename+'.flux'
00165                 #       pixels=result+'/'+fluxscale_image
00166                 #       ia.calc(pixels=pixels)
00167                 #       ia.close()
00168 
00169 
00170         except Exception, instance:
00171                 print '*** Error *** ',instance
00172                 raise Exception, instance
00173