set user     "AST"
set source   "Abell1205"
set srclower "abell1205"
set date     "20110415"
set cycle    "8"
set mode     "Off"
set toppath  "/sps/baoradio/AmasNancay/${user}"
set path     "${toppath}/${source}/dab"
set ppf      "0"

# Read DAB calibration in Jansky (Voie E / PKS1127-14)
#-----------------------------------------------------
n/read nt_dabjy ${toppath}/DAB_Jy.dat " " dabfreq,0 dabjy,1
newwin 1 1
n/pl nt_dabjy.dabjy%dabfreq 1 ! "marker=star,4 magenta notit"
settitle "DAB Voie E / PKS1127-14. Juin 2011"


# Read DAB calibration in A.U. (from BAORadio data)
#--------------------------------------------------
openppf ${path}/dabchrom_${date}_${srclower}_${mode}.ppf
ls

## Baseline
del b0 b1 b0f b1f
objaoper /home/base row 0 b0
objaoper /home/base row 1 b1
objaoper /home/basefilt row 0 b0f
objaoper /home/basefilt row 1 b1f

## First DAB
del dab01 dab11 dab01f dab11f dab01n dab11n
objaoper /home/firstdab row 0 dab01
objaoper /home/firstdab row 1 dab11
objaoper /home/firstdabfilt row 0 dab01f
objaoper /home/firstdabfilt row 1 dab11f
objaoper /home/firstdelta row 0 dab01d
objaoper /home/firstdelta row 1 dab11d

## Second DAB
del dab02 dab12 dab02f dab12f dab02n dab12n
objaoper /home/seconddab row 0 dab02
objaoper /home/seconddab row 1 dab12
objaoper /home/seconddabfilt row 0 dab02f
objaoper /home/seconddabfilt row 1 dab12f
objaoper /home/seconddelta row 0 dab02d
objaoper /home/seconddelta row 1 dab12d

## Mean(delta(first,second DAB)) per channel
c++exec TVector<r_4> meandab0d(dab01.NElts()); KeepObj(meandab0d); \
        TVector<r_4> meandab1d(dab01.NElts()); KeepObj(meandab1d); \
        meandab0d=dab01d+dab02d; meandab0d/=2.; \
        meandab1d=dab11d+dab12d; meandab1d/=2.;

## Save Mean(delta(first,second DAB)) per channel into a ppf file
#if ( ppf == "1" ) then
#  del daboff
#  c++exec TMatrix<r_4> daboff(2,dab01.NElts()); KeepObj(daboff); \
#          daboff.Row(0)=meandab0d; daboff.Row(1)=meandab1d;
#  saveppf daboff "${toppath}/${source}/daboff_${date}_${srclower}_calibcycle8.fits"
#endif 

# Plot results
#--------------
graphicatt "xylimits=1250,1500,0,3"

## Baseline
newwin 1 1
plot2d b0 (n/8192)*250+1250 val n>0 "blue cpts nstat notit"
plot2d b1 (n/8192)*250+1250 val n>0 "same red cpts nstat notit"
settitle "Baseline raw $source $date $mode Ch0 blue Ch1 red"

newwin 1 1
plot2d b0f (n/8192)*250+1250 val n>0 "blue cpts nstat notit"
plot2d b1f (n/8192)*250+1250 val n>0 "same red cpts nstat notit"
settitle "Baseline filt $source $date $mode Ch0 blue Ch1 red"

## First DAB
#newwin 1 1
#plot2d dab01 (n/8192)*250+1250 val n>0 "blue cpts nstat notit"
#plot2d dab11 (n/8192)*250+1250 val n>0 "same red cpts nstat notit"
#settitle "First DAB raw $source $date $mode Ch0 blue Ch1 red"

newwin 1 1
plot2d dab01f (n/8192)*250+1250 val n>0 "blue cpts nstat notit"
plot2d dab11f (n/8192)*250+1250 val n>0 "same red cpts nstat notit"
settitle "First DAB filt $source $date $mode Ch0 blue Ch1 red"

## Second DAB
#newwin 1 1
#plot2d dab02 (n/8192)*250+1250 val n>0 "blue cpts nstat notit"
#plot2d dab12 (n/8192)*250+1250 val n>0 "same red cpts nstat notit"
#settitle "Second DAB raw $source $date $mode Ch0 blue Ch1 red"

newwin 1 1
plot2d dab02f (n/8192)*250+1250 val n>0 "blue cpts nstat notit"
plot2d dab12f (n/8192)*250+1250 val n>0 "same red cpts nstat notit"
settitle "Second DAB filt $source $date $mode Ch0 blue Ch1 red"

## Delta=DAB-baseline + Calibrated DAB in Jansky
graphicatt "xylimits=1100,1600,0,30"
newwin 1 1
n/pl nt_dabjy.dabjy%dabfreq 1 ! "marker=star,4 turquoise notit nstat"
plot2d dab01d (n/8192)*250+1250 val*18 n>0 "same blue cpts nstat notit"
plot2d dab11d (n/8192)*250+1250 val*22 n>0 "same red cpts nstat notit"
settitle "First DAB delta $source $date $mode"
setaxelabels "Freq (MHz)" "Intensity" "bigboldfont"
addline 1420 4 1430 4 "turquoise line=solid,3"
addline 1420 3 1430 3 "blue line=solid,3"
addline 1420 2 1430 2 "red line=solid,3"
addtext 1440 4 "DAB Voie E (Jy)" "bigfont"
addtext 1440 3 "Ch0*18 (A.U.)" "bigfont"
addtext 1440 2 "Ch1*22 (A.U.)" "bigfont"

newwin 1 1
n/pl nt_dabjy.dabjy%dabfreq 1 ! "marker=star,4 turquoise notit nstat"
plot2d dab02d (n/8192)*250+1250 val*18 n>0 "same blue cpts nstat notit"
plot2d dab12d (n/8192)*250+1250 val*22 n>0 "same red cpts nstat notit"
settitle "Second DAB delta $source $date $mode"
setaxelabels "Freq (MHz)" "Intensity" "bigboldfont"
addline 1420 4 1430 4 "turquoise line=solid,3"
addline 1420 3 1430 3 "blue line=solid,3"
addline 1420 2 1430 2 "red line=solid,3"
addtext 1440 4 "DAB Voie E (Jy)" "bigfont"
addtext 1440 3 "Ch0*18 (A.U.)" "bigfont"
addtext 1440 2 "Ch1*22 (A.U.)" "bigfont"

## Delta=Mean(DAB-baseline) + Calibrated DAB in Jansky
graphicatt "xylimits=1100,1600,0,30"
newwin 1 1
n/pl nt_dabjy.dabjy%dabfreq 1 ! "marker=star,4 turquoise notit nstat"
plot2d meandab0d (n/8192)*250+1250 val*18 n>0 "same blue cpts nstat notit"
plot2d meandab1d (n/8192)*250+1250 val*23 n>0 "same red cpts nstat notit"
settitle "Mean DAB delta $source $date $mode"
setaxelabels "Freq (MHz)" "Intensity" "bigboldfont"
addline 1420 4 1430 4 "turquoise line=solid,3"
addline 1420 3 1430 3 "blue line=solid,3"
addline 1420 2 1430 2 "red line=solid,3"
addtext 1440 4 "DAB Voie E (Jy)" "bigfont"
addtext 1440 3 "Ch0*18 (A.U.)" "bigfont"
addtext 1440 2 "Ch1*23 (A.U.)" "bigfont"