# split out the fields from a mosaic

include 'ms.g'
myms:=ms('orion_source.ms',readonly=F);

for (k in 1:9){
msname:=spaste('FIELD_',as_string(k),'.ms');
myms.split(outputms=msname,fieldids=[k],spwids=[1,2], whichcol='CORRECTED_DATA');
}

myms.done();

# now image all of the fields (in separate MSs) together

include 'imager.g'
include 'calibrater.g'
include 'logger.g'
im:=imager();

for (k in 1:9){
msname:=spaste('FIELD_',as_string(k),'.ms');
im.setdata(mode="none" ,  spwid=[1, 2] ,
        fieldid=1, msname=msname);
}

# set cell and image size, reference fieldid = 5, mosaic center:
ok:=im.setimage(nx=300, ny=300, cellx='2arcsec', celly='2arcsec',
        stokes="I", mode="mfs",spwid=[1],fieldid=1,doshift=T,
        phasecenter=[m0=[value=83.8227917, unit="deg" ],m1=[value=-5.38521944, unit="deg" ], 
	type="direction",refer="J2000"]);
# set the voltage pattern for the VLA antennas:
ok:=im.setvp(dovp=T);
# use robust weighting:
ok:=im.setoptions(ftmachine='mosaic');
ok:=im.setmfcontrol(stoplargenegatives=-1, cyclefactor=3, cyclespeedup=500)
ok:=im.weight(type="briggs", rmode="norm", robust=-1,mosaic=T);
im.make('orion.B2');
  # Weighting MS: IMAGING_WEIGHT column will be changed
  # Briggs weighting: sidelobes will be suppressed over full image
  #                  : using 300 pixels in the uv plane
# NOTE: use defaults on cyclefactor and cyclespeedup so no need to
#  run the setup function myimager.setmfcontrol
#
# Run MEM with full image allowed in clean box:
# view Orion mask
#
ok:=im.mem(algorithm="mfentropy" , niter=50,
        sigma='4mJy', targetflux='180Jy',constrainflux=F,
        model="orion.B2", mask="orion.mask6" , displayprogress=T,
        image="orion.B2.restored" , residual="orion.B2.residual" );