#
#
include 'ms.g';
include 'autoflag.g';
include 'calibrater.g';
include 'imager.g'
#
mset:=fitstoms(msfile='n1333.ms',
fitsfile='/home/mzwaan/cdrom/Data/N1333/N1333_1.UVFITS');
mset.summary();
mset.done();
#
af:=autoflag('n1333.ms');
af.setdata();
af.setselect(autocorr=T);
af.run();
af.setselect(ant=[22],spwid=[2]);
af.run();
# bandpasses for ant 9 look bad
af.setselect(ant=[9],spwid=[1,2]);
af.run();
af.setselect(field=[2:6,8:12],clip=[expr="+ ABS RR LL",max=.75]);
af.run();
af.setselect(timerng="2003/05/02/21:44 2003/05/02/21:50",spwid=[1]);
af.run();
af.setselect(timerng="2003/05/02/21:40 2003/05/02/21:51",spwid=[1]);
af.run();
af.setselect(timerng="2003/05/02/21:56 2003/05/02/21:58",spwid=[2]);
af.run();
af.setselect(timerng="2003/05/02/21:55 2003/05/02/22:20",spwid=[2]);
af.run();
af.setselect(chan=[1,3,61,63],spwid=[1,2]);
af.run();
af.setselect(timerng="2003/05/02/00:00 2003/05/02/19:30",ant=[14]);
af.run();
af.setselect(timerng="2003/05/02/00:00 2003/05/02/19:20",ant=[15]);
af.run();
af.setselect(timerng="2003/05/02/00:00 2003/05/02/20:00",ant=[18]);
af.run();
af.setselect(timerng="2003/05/02/00:00 2003/05/02/19:00",ant=[22]);
af.run();
af.setselect(timerng="2003/05/03/00:40 2003/05/03/02:00",ant=[4]);
af.run();
af.done;
#
#
# SETJY
#
imgr:=imager(filename='n1333.ms');
imgr.setjy(fieldid=14, spwid=2);
imgr.setjy(fieldid=13, spwid=1);
imgr.done();
#
# SOLVE FOR GAINS
# FIRST DO THE BANDPASS CALIBRATORS:
#
cal:=calibrater(filename='n1333.ms');
cal.reset();
cal.setdata(msselect='FIELD_ID in [15:16]',mode='channel',start=3,nchan=55,
uvrange=[0,50]);
cal.setapply(type='GAINCURVE',t=-1);
cal.setsolve(type='G',t=0,refant=26,table='n1333.gcal',append=F);
cal.solve();
#
# FIND BANDPASSES
#
cal.reset();
cal.setdata(msselect='FIELD_ID in [15:16]',uvrange=[0,50]);
cal.setapply(type='GAINCURVE',t=-1);
cal.setapply(type='G',t=0.0,table='n1333.gcal');
cal.setsolve(type='B',t=0,refant=26,table='n1333.bcal',append=F);
cal.state();
cal.solve();
#
# PLOT BANDPASSES
#
cal.plotcal(tablename='n1333.bcal',plottype='AMP',multiplot=T,nx=3,ny=3);
#
#
# FIND BANDPASSES, NOW WITH BPOLY
#
#cal.reset();
#cal.setdata(msselect='FIELD_ID in [15]');
#cal.setapply(type='GAINCURVE',t=-1);
#cal.setapply(type='G',t=0.0,table='n1333.gcal');
#cal.setsolvebandpoly(table='n1333.bpoly1',degamp=5,degphase=7,append=F);
#cal.state();
#cal.solve();
#cal.reset();
#cal.setdata(msselect='FIELD_ID in [16]');
#cal.setapply(type='GAINCURVE',t=-1);
#cal.setapply(type='G',t=0.0,table='n1333.gcal');
#cal.setsolvebandpoly(table='n1333.bpoly2',degamp=5,degphase=7,append=F);
#cal.state();
#cal.solve();
#
#shell('gv n1333.bpoly1.ps &');
#shell('gv n1333.bpoly2.ps &');
#
#
# FIND GAINS FOR SECONDARY CALS
# APPLY BANDPASSES
#
cal.reset();
cal.setdata(msselect='FIELD_ID IN [1,7]',mode='channel',start=3,nchan=55);
cal.setapply(type='GAINCURVE',t=-1);
cal.setapply(type='B',t=0.0,table='n1333.bcal');
cal.setsolve(type='G',t=0,refant=26,table='n1333.gcal',append=T);
cal.solve();
#
# PLOT PHASES OF SECONDARIES
#
cal.plotcal(tablename='n1333.gcal',plottype='PHASE',fields=[1,7],
nx=3,ny=3,multiplot=T);
cal.plotcal(tablename='n1333.gcal',plottype='AMP',fields=[1,7],
nx=3,ny=3,multiplot=T);
#
# find gains for secondary cals, but now do a spline fit
# apply bandpasses, can't apply G, not implemented
# does not work because antenna 13 has no data and it stops
# after baseline 1-13. The setdata does not help.
#
#cal.reset();
#cal.setdata(msselect=
#'FIELD_ID IN [1,7] && ANTENNA1 in [1:12,14:27] && ANTENNA2 in [1:12,14:27]',
#mode='channel',start=3,nchan=55);
#cal.setapply(type='GAINCURVE',t=-1);
#cal.setapply(type='B',t=0.0,table='n1333.bcal');
#cal.setsolvegainspline(mode='PHAS',splinetime=10000,refant=26,
#table='n1333.cal_ph.gspl',append=F);
#cal.solve();
#
#
# FIND GAINS FOR FLUX CALS
#
cal.reset();
cal.setdata(msselect='FIELD_ID IN [13,14]',mode='channel',start=3,nchan=55);
cal.setapply(type='GAINCURVE',t=-1);
cal.setapply(type='B',t=0.0,table='n1333.bcal');
cal.setsolve(type='G',t=0,refant=26,table='n1333.gcal',append=T);
cal.solve();
#
# TRANSFER FLUXES TO SECONDARIES
#
cal.fluxscale(tablein='n1333.gcal',tableout='n1333.fluxcal',
reference="0542+498_2 0542+498_1",transfer="0336+323_2 0336+323_1");
#
#Flux density for 0336+323_1 in SpW=1 is: 2.24716 +/- 0.00662403 (SNR = 339.244
#Flux density for 0336+323_2 in SpW=2 is: 2.38358 +/- 0.00814374 (SNR = 292.689
#
# FIRST APPLY CALIBRATION TO THE CALIBRATOR:
#
cal:=calibrater(filename='n1333.ms');
cal.reset();
cal.setdata(msselect='FIELD_ID IN [1,7] && SPECTRAL_WINDOW_ID in [1,2]');
cal.setapply(type='B',t=0.0,table='n1333.bcal',select='');
cal.setapply(type='GAINCURVE',t=-1);
cal.setapply(type='G',t=0.0,table='n1333.fluxcal',select='FIELD_ID in [1,7]');
cal.correct();
#
# THEN TO THE FIELDS:
#
cal.reset();
cal.setdata(msselect='FIELD_ID IN [2:6,8:12] && SPECTRAL_WINDOW_ID in [1,2]');
cal.setapply(type='GAINCURVE',t=-1);
cal.setapply(type='B',t=0.0,table='n1333.bcal',select='');
cal.setapply(type='G',t=0.0,table='n1333.fluxcal',select='FIELD_ID in [1,7]');
cal.state();
cal.correct();
#
cal.done();
#
#
# DONE WITH CALIBRATION
#
#
# SPLIT THE DATA
#
myms:= ms("n1333.ms",readonly=F);
myms.split('n1333cal.ms',fieldids=[2:6,8:12],spwids=[1,2],
whichcol='CORRECTED_DATA',nchan=[63,63],start=[1,1],step=[1,1]);
myms.done();
#
af:=autoflag('n1333cal.ms');
af.setdata();
af.setselect(field=[2:6,8:12],clip=[expr="+ ABS RR LL",max=7]);
af.run();
af.done;
#
#
# IMAGE
#
imgr:=imager('n1333cal.ms');
#
imgr.setdata(mode='channel',fieldid=[1:10],nchan=[63,63],start=[1,1],
step=[1,1],spwid=[1,2]);
imgr.setimage(nx=640,ny=256,cellx='0.5arcsec',celly='0.5arcsec',stokes='I',
mode='channel',nchan=18,start=3,step=5,fieldid=[2],spwid=[1,2]);
imgr.weight(type='briggs',rmode='norm',robust=0.5,mosaic=T);
imgr.setvp(dovp=T);
imgr.setoptions(ftmachine='mosaic',padding=1.0);
imgr.setmfcontrol(scaletype='SAULT',minpb=0.07,constpb=0.4);
imgr.clean(algorithm='mfclark',niter=1000,gain=0.15,model='model_mos',
residual='residual_mos',image='restored_mos',interactive=F);
#
#