#! /bin/csh -f
##############################################
#
#2016-10-24: JHZ
#2016-10-28: JHZ
#2016-11-01: JHZ
#2016-11-09: JHZ
#2016-11-15: JHZ
#2016-11-18: JHZ
#2016-12-27: JHZ
#2017-01-13: JHZ
#2017-07-06: JHZ
#version:
#miriad-sma3.3.4
#miriad-sma4.0.0
#miriad-sma4.0.1
#miriad-sma4.0.2
#miriad-sma4.0.3
#miriad-sma5.0.0
#miriad-sma5.0.1
#miriad-sma5.0.2
#miriad-sma5.0.3 and higher versions
#
#test CSDI for SMA polarization data
#using the stokes i, q, u, v produced from SMA
#and deep cleaned with the newly implemented CSDI
#the calibrated polarization uv data imported
#from case 12 of the examples given in
#SMA Miriad scripts
#
#This script provides a comparison between
#outputs of from CSDI and clean.
#
#Also a demo for displaying polarization vecotors, intensity
# along with Stokes Q, U, V, I
########## define variables ##########
#
set source = 3c273
set sb = usb
set polvis = pp3c273_flx.usb.blc
set nw = 48
set imsize = 512
set cell = 0.1
set disreg='arcsec,box(-8,-8,8,8)'
set clreg='arcsec,box(-3,-3,3,3)'
set pdisreg='arcsec,box(-3.5,-3.5,3.5,3.5)'
set statreg='arcsec,box(-10,-10,-8,-8)'
#
########## process handles ##########
#
#goto INVERT # FFT
#goto CLEAN # classical clean
#goto CSDI # steering clean
#goto RESTOR # convolv to a clean beam
#goto POLVEC # display map of polarization vectors
#goto DISPIQUV # display four stokes contours
goto DISPVEC # display pol vector, intensity along with stokes Q,U,V
#
########## process with Miriad routines ##########
INVERT:
\rm -fr $source.$sb.i $source.$sb.q $source.$sb.u $source.$sb.v $source.$sb.b
invert vis=$polvis beam=$source.$sb.b \
imsize=$imsize cell=$cell \
map=$source.$sb.i,$source.$sb.q,$source.$sb.u,$source.$sb.v \
options=systemp sup=0 \
stokes=i,q,u,v robust=2 line=chan,1,1,$nw,1
$<
cgdisp in=$source.$sb.i,$source.$sb.i type=p,c \
device=/xw \
labtyp=arcsec nxy=1,1 slev=p,10 \
levs1=2,3,4,5,6,7,8,9 \
options=wedge,mirror region=$disreg
exit
CLEAN:
\rm -fr $source.$sb.icln $source.$sb.qcln $source.$sb.ucln $source.$sb.vcln
\rm -fr $source.$sb.icmp $source.$sb.qcmp $source.$sb.ucmp $source.$sb.vcmp
echo 'classic clean - i'
clean map=$source.$sb.i beam=$source.$sb.b out=$source.$sb.icmp \
region=$clreg gain=0.05 niters=150
$<
restor map=$source.$sb.i beam=$source.$sb.b out=$source.$sb.icln \
model=$source.$sb.icmp
$<
cgdisp in=$source.$sb.icln,$source.$sb.icln type=p,c \
device=/xw \
labtyp=arcsec nxy=1,1 slev=p,10 \
levs1=2,3,4,5,6,7,8,9 \
options=wedge,mirror region=$disreg
$<
imstat in=$source.$sb.icln region=$statreg
$<
echo 'classic clean - q'
clean map=$source.$sb.q beam=$source.$sb.b out=$source.$sb.qcmp \
region=$clreg gain=0.05 niters=150
$<
restor map=$source.$sb.q beam=$source.$sb.b out=$source.$sb.qcln \
model=$source.$sb.qcmp
$<
cgdisp in=$source.$sb.qcln,$source.$sb.icln type=p,c \
device=/xw \
labtyp=arcsec nxy=1,1 slev=p,10 \
levs1=2,3,4,5,6,7,8,9 \
options=wedge,mirror region=$disreg
$<
imstat in=$source.$sb.qcln region=$statreg
$<
echo 'classic clean - u'
clean map=$source.$sb.u beam=$source.$sb.b out=$source.$sb.ucmp \
region=$clreg gain=0.05 niters=150
$<
restor map=$source.$sb.u beam=$source.$sb.b out=$source.$sb.ucln \
model=$source.$sb.ucmp
$<
cgdisp in=$source.$sb.ucln,$source.$sb.icln type=p,c \
device=/xw \
labtyp=arcsec nxy=1,1 slev=p,10 \
levs1=2,3,4,5,6,7,8,9 \
options=wedge,mirror region=$disreg
$<
imstat in=$source.$sb.ucln region=$statreg
$<
echo 'classic clean - v'
clean map=$source.$sb.v beam=$source.$sb.b out=$source.$sb.vcmp \
region=$clreg gain=0.05 niters=150
$<
restor map=$source.$sb.v beam=$source.$sb.b out=$source.$sb.vcln \
model=$source.$sb.vcmp
$<
cgdisp in=$source.$sb.vcln,$source.$sb.icln type=p,c \
device=/xw \
labtyp=arcsec nxy=1,1 slev=p,10 \
levs1=2,3,4,5,6,7,8,9 \
options=wedge,mirror lines=1,5,5 region=$disreg
$<
imstat in=$source.$sb.vcln region=$statreg
exit
DISPIQUV:
cgdisp in=$source.$sb.vcln,$source.$sb.qcln,$source.$sb.ucln,$source.$sb.icln \
type=c,c,c,p \
device=/xw \
labtyp=arcsec nxy=1,1 slev=p,10 \
levs1=2,3,4,5,6,7,8,9 \
options=wedge,mirror lines=1,2,4,8 region=$disreg
exit
$<
CSDI:
echo 'Steer deep clean -'
\rm -rf ${source}.${sb}.q.csdi ${source}.${sb}.u.csdi
csdi map=${source}.${sb}.q,${source}.${sb}.u \
beam=3c273.usb.b \
out=${source}.${sb}.q.csdi,${source}.${sb}.u.csdi \
gain=0.05 niters=500
echo 'done CSDI'
RESTOR:
\rm -rf ${source}.${sb}.csdi.qcln ${source}.${sb}.csdi.ucln
$<
restor map=$source.$sb.q beam=$source.$sb.b out=$source.$sb.csdi.qcln \
model=$source.$sb.q.csdi
$<
restor map=$source.$sb.u beam=$source.$sb.b out=$source.$sb.csdi.ucln \
model=$source.$sb.u.csdi
$<
echo 'CSDI - q'
cgdisp in=$source.$sb.csdi.qcln,$source.$sb.icln type=p,c \
device=/xw \
labtyp=arcsec nxy=1,1 slev=p,10 \
levs1=2,3,4,5,6,7,8,9 \
options=wedge,mirror region=$disreg
$<
imstat in=$source.$sb.csdi.qcln region=$statreg
$<
echo 'CSDI - u'
cgdisp in=$source.$sb.csdi.ucln,$source.$sb.icln type=p,c \
device=/xw \
labtyp=arcsec nxy=1,1 slev=p,10 \
levs1=2,3,4,5,6,7,8,9 \
options=wedge,mirror region=$disreg
$<
imstat in=$source.$sb.csdi.ucln region=$statreg
exit
$<
POLVEC:
set irms=0.01
set qurms=0.005
\rm -fr $source.$sb.pi $source.$sb.pie
\rm -fr $source.$sb.pm $source.$sb.pme
\rm -fr $source.$sb.pa $source.$sb.pae
$<
impol in=$source.$sb.csdi.qcln,$source.$sb.csdi.ucln,$source.$sb.icln \
poli=$source.$sb.pi,$source.$sb.pie \
polm=$source.$sb.pm,$source.$sb.pme \
pa=$source.$sb.pa,$source.$sb.pae \
sigma=$qurms,$irms sncut=3,10
$<
echo 'polarization vectors -'
cgdisp in=$source.$sb.pi,$source.$sb.pm,$source.$sb.pa \
type=pix,amp,ang device=/xs \
labtyp=arcsec options=mirror,blacklab \
vecfac=2.0,3.0 cols1=0.5 \
range=0,0,lin,1 region=$pdisreg
$<
echo 'error of polarization vectors -'
cgdisp in=$source.$sb.pie,$source.$sb.pme,$source.$sb.pae \
type=pix,amp,ang device=/xs \
labtyp=arcsec options=mirror,blacklab \
vecfac=2.0,3.0 cols1=0.5 \
range=0,0,lin,1 region=$pdisreg
$<
echo 'in comparison to classical method -'
\rm -fr $source.$sb.pi $source.$sb.pie
\rm -fr $source.$sb.pm $source.$sb.pme
\rm -fr $source.$sb.pa $source.$sb.pae
$<
impol in=$source.$sb.qcln,$source.$sb.ucln,$source.$sb.icln \
poli=$source.$sb.pi,$source.$sb.pie \
polm=$source.$sb.pm,$source.$sb.pme \
pa=$source.$sb.pa,$source.$sb.pae \
sigma=$qurms,$irms sncut=3,10
$<
echo 'polarization vectors -'
cgdisp in=$source.$sb.pi,$source.$sb.pm,$source.$sb.pa \
type=pix,amp,ang device=/xs \
labtyp=arcsec options=mirror,blacklab \
vecfac=2.0,3.0 cols1=0.5 \
range=0,0,lin,1 region=$pdisreg
$<
echo 'error of polarization vectors -'
cgdisp in=$source.$sb.pie,$source.$sb.pme,$source.$sb.pae \
type=pix,amp,ang device=/xs \
labtyp=arcsec options=mirror,blacklab \
vecfac=2.0,3.0 cols1=0.5 \
range=0,0,lin,1 lines=1,9 region=$pdisreg
exit
DISPVEC:
cgdisp in=$source.$sb.pi,$source.$sb.pm,$source.$sb.pa,$source.$sb.ucln,$source.$sb.qcln,$source.$sb.csdi.vcln \
type=pix,amp,ang,c,c,c device=/xs \
labtyp=arcsec options=mirror,blacklab \
vecfac=1.5,4.0 cols1=3.0 cols2=1.0 cols3=7.0 \
range=0,0,lin,1 region=$pdisreg \
lines=1,7,5,7,9 \
slev=p,10,p,10,p,10 \
levs1=4,8,16,28,44 \
levs2=4,5,6,7,8,9 \
levs3=4,5,6,7,8,9
###########################
#lines(1) = 1 frame & lable
#lines(2) = 7 stokes u => cols1=3.0 green
#lines(3) = 5 stokes q => cols2=1.0 white
#lines(5) = 7 stokes v => cols3=7.0 yellow
#lines(7) = 9 pol vector => orange
###########################
terminate:
echo 'end'
Figure from the process of DISPVEC