#! /bin/csh -f










# Some shortcuts
source globals.csh

# Define defaults
set antpos=antpos.24sep09       # Antenna positions file
set bpref=1                     # Antenna used as refant for BP
set cell=5
set imsize=256
set rawDir=DATA/

### Usage
if ($#argv<2) then
    echo "Usage: "
    echo "  $0 <date> <tracks> <start> <express?> <bpcal>"
    echo "    <date> is date string that will be used to name tracks (e.g. sep24)"
    echo "    <tracks> is double-quoted, comma-separated list of miriad files"
    echo "    <start> is a startpoint in the script"
    echo "      GET, WIN, BASELINE, FLAG, LINECAL, BP"
    echo "    <express> is y/n. Default: n. If 'y', skip plotting steps"
    echo "    <bpcal> is the source to use for bandpass calibration"
    exit 1
endif

set date=$1
set tracks=$2   # List of tracks, quoted, time-ordered
if ($#argv>2) then
    set start=$3
    set exp='n'
    if ($#argv>3) then
	if ($4 == 'y') set exp='y'
    endif
    if ($#argv>4) then
	set bp=$5
	set dobp=1
    else
	set dobp=0
    endif
else
    set start='GET'
endif


### Prepare log file
set log=$date'_CAL.log'
if ($start == 'GET') then
#    echo here; set ans=$<
    ifrm $log
else
    set logline=`grep -n "###>$start" $log | cut -f1 -d':'`
    @ logline--
    mv $log $log.2 
    head -$logline $log.2 > $log
    rm $log.2
endif
set flg=$date.flags

### Summarize input
Tout "#################################################################" $log
Tout "##### Inputs: $argv " $log
Tout "#####   Date: $date " $log
Tout "#####   Tracks: $tracks " $log
Tout "#####          " $log
Tout "#################################################################" $log

goto $start

##################################################
GET:
echo "###>GET" >> $log
### Read in data file(s) and concatenate if multiple
ifrm $date.mir
set thisDir=`pwd`
cd $rawDir
Tout "uvcat vis=$tracks out=$thisDir/$date.mir" $thisDir/$log
uvcat vis=$tracks out=$thisDir/$date.mir
cd $thisDir
### Grab start and end of track
uvindex vis=$date.mir log=$date.index
set startStr = \
  `\grep -A3 "Antennas Spectral" $date.index | tail -1 | cut -f1 -d' '`
set endStr = \
  `\grep -B1 "Total number of records" $date.index | head -1 | cut -f1 -d' '`
Tout "##### Track Start: $startStr" $log
Tout "##### Track End:   $endStr" $log
./parseIndex.csh $date
Tout "##### Sources: " $log
cat $date.srcs
cat $date.srcs >> $log

##################################################
### Check that systemp is different between sidebands
# If equal, sideband ratio was wrong. Fix Tsys
# This could be implemented as a search for 0.0000000 and 10000.000

### Check that PAM atten is what it's supposed to be

##################################################
WIN:
echo "###>WIN" >> $log
### Divide by autocorrelations
ifrm $date.auto
Tout "##### Correcting by autocorrelations" $log
uvcal vis=$date.mir out=$date.auto options=fxcal
### write out windows 1-3, drop autos
ifrm $date.lsb
uvcat vis=$date.auto out=$date.lsb select='window(1,2,3),-auto,-source(noise)'

##################################################
BASELINE:
echo "###>BASELINE" >> $log
### Test whether correction is needed
set startStr = \
  `\grep -A3 "Antennas Spectral" $date.index | tail -1 | cut -f1 -d' '`
# convert day to number (days from Sep 0th)
set dayNum=`echo $startStr | cut -c3-7 |sed -e's/SEP/0  /' -e's/OCT/30 /' |awk '{print $1+$2}'`
# convert hours to fractional day
set hrNum=`echo $startStr | cut -c9-16 |sed -e's/:/ /' |awk '{print ($1+($2+$3/60.)/60.)/24.}'`
# Now compare to 09SEP25:09:50:15, when baselines were changed
set before=`echo $dayNum $hrNum 25.409722 | awk '{print ($3-$2-$1)>=0}'`
if ($before == '1') then
    ### Correct baselines to best available knowledge
    ifrm $date.bas0
    Tout "##### Updating antenna positions using file: $antpos" $log
    uvedit vis=$date.lsb out=$date.bas0 apfile=$antpos
    cat $antpos >> $log
    echo "############################################################" >>$log
else
    Tout "##### No baseline correction required for date $startStr" $log
    rmcp $date.lsb $date.bas
endif

### Now check whether NIA correction is required
# NIA terms go to correct values: 09SEP25:09:50:15  (as above)
# Before that, went to correct values with wrong sign 09SEP21:19:32:22
# Before that, no NIA in RTS
if ($before == '1') then
    set dblNIAtime=21.814144
    set before2=`echo $dayNum $hrNum $dblNIAtime |awk '{print ($3-$2-$1)>=0}'`
    if ($before2 == '1') then
	### Here we need to add single NIA correction
	Tout "##### Date $startStr is before 09SEP21:19:32:22. No NIA in RTS" $log
	Tout "##### Correcting NIA" $log
	ifrm $date.nia
	mv $date.bas0 $date.nia
	axcor vis=$date.nia out=$date.bas0 \
	    axoff=0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-0.0109,0.0035,-0.0070,-0.0011,-0.0118,0.0098,0.0041,-0.0099
	rm -r $date.nia
    else
	### Here we need to add double NIA correction (sign wrong)
	Tout "##### Date $startStr is before 09SEP25:09:50:15 and after 09SEP21:19:32:22. NIA incorrect in RTS" $log
	Tout "##### Correcting NIA" $log
	ifrm $date.nia
	mv $date.bas0 $date.nia
	axcor vis=$date.nia out=$date.bas0 \
	    axoff=0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-0.0218,0.0070,-0.0140,-0.0022,-0.0236,0.0196,0.0082,-0.0198
	rm -r $date.nia
    endif
else
    Tout "##### No NIA correction required for date $startStr" $log
endif

### Regenerate wideband
if ($before == '1') then
    Tout "##### Regenerating the wideband data" $log
    ifrm $date.bas
    uvcal vis=$date.bas0 out=$date.bas options=avechan
    rm -r $date.bas0
endif


##################################################
FLAG:
echo "###>FLAG" >> $log
rmcp $date.bas $date.flg
### flag Tsys jumps
## For now we just have to find these by hand and flag them

### Check PAM attenuations
./checkPamAttens.sh $tracks | tee -a $log
if !($exp == 'y') then
    echo "   Press enter to continue" ; set ans=$<
endif

if !($exp == 'y') then
    echo ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"
    echo ">>>>> Plotting Tsys data"
    echo ">>>>> Find bad points, insert in flag file $flg"
    echo ">>>>>   After completion of this step, $flg will be sourced"
    echo "   Press enter to continue" ; set ans=$<
    ./flagPlot.csh $date $startStr $endStr
    echo ">>>>> Press enter when $flg is up to date" ; set ans=$<
endif

FLAGCONT:
echo "###>FLAGCONT" >> $log
Tout "##### Importing flags from $flg" $log
cat $flg >> $log
echo "############################################################" >> $log
source $flg

### csflag
csflag vis=$date.flg sza=true carma=false

##################################################
LINECAL:
echo "###>LINECAL" >> $log
### Linelength Calibration
ifrm $date.lcg
./linecalx vis=$date.flg out=$date.lcg

### Display LL phase
if !($exp == 'y') then
    echo ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"
    echo ">>>>> Displaying linelength phase "
    echo "   Press enter to continue" ; set ans=$<
    gpplt vis=$date.lcg yaxis=phase yrange=-180,180 options=wrap,dots \
    device=/xs nxy=3,5
    set ans=$<
endif

### Apply LL phase
Tout "##### Applying linelength calibration " $log
gpcopy vis=$date.lcg out=$date.flg

### Create calibrated file
ifrm $date.lc
uvaver vis=$date.flg out=$date.lc

### Display correction
if !($exp == 'y') then
    echo ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"
    echo ">>>>> Displaying linelength correction"
    echo "  >>> white:  before correction"
    echo "  >>> yellow: after correction"
    echo "   Press enter to continue" ; set ans=$<
    uvplt vis=$date.flg,$date.lc axis=ti,ph nxy=5,5 device=/xw options=nocal \
    select="-auto,-source($clusters,noise)"
    echo "   Press enter to continue" ; set ans=$<
endif

##################################################
BP:
echo "###>BP" >> $log
### Bandpass Calibration
Tout "##### Bandpass Calibration" $log

### ensure that bp source was selected
if !($?bp) then
    echo ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"
    echo ">>>>> Plotting BP calibrator candidates."
    echo "   Press enter to continue" ; set ans=$<
    smauvplt vis=$date.lc device=/xw options=nob select="-source(jupiter,$clusters)"
    echo -n " Select bandpass calibrator: " ; set bp=$<
endif
Tout "### Source $bp" $log
ifrm $date.bp
uvcat vis=$date.lc out=$date.bp select="source($bp)"
mfcal vis=$date.bp interval=4 refant=$bpref
if !($exp == 'y') then
    gpplt vis=$date.bp device=/xw nxy=5,3 yaxis=amp
    echo ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"
    echo ">>>>> Gain amps for BP calibrator $bp"
    echo "   Press enter to continue" ; set ans=$<
    gpplt vis=$date.bp device=/xw nxy=5,3 yaxis=pha
    echo ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"
    echo ">>>>> Gain phase for BP calibrator $bp"
    echo "   Press enter to continue" ; set ans=$<
    gpplt vis=$date.bp device=/xw nxy=5,3 options=bandpass yrange=0,2
    echo ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"
    echo ">>>>> Bandpass amp for BP calibrator $bp"
    echo "   Press enter to continue" ; set ans=$<
    gpplt vis=$date.bp device=/xw nxy=5,3 options=bandpass \
    yaxis=pha yrange=-180,180
    echo ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"
    echo ">>>>> Bandpass phase for BP calibrator $bp"
    echo "   Press enter to continue" ; set ans=$<
endif
## Apply BP, regen wide, scale cal flux to unity for display purposes
ifrm $date.bpPL
uvcal vis=$date.bp out=$date.bpPL options=avechan
mfcal vis=$date.bpPL options=nopass flux=1,34,0
if !($exp == 'y') then
    echo ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"
    echo ">>>>> Calibrator BP after calibration (should be very flat)"
    echo "   Press enter to continue" ; set ans=$<
    smauvspec vis=$date.bpPL device=/xw nxy=5,5 axis=freq,both \
    interval=1000 yrange=0,2
    echo "   Press enter to continue" ; set ans=$<
endif

BPTEST:
### Now display another source as a test
if !($exp == 'y') then
    smauvplt vis=$date.lc device=/xw line=ch,1,1,15 options=nob select='-source(mars,jupiter)'
    echo ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"
    echo ">>>>> Select target to test BP solution"
    echo ">>>>> Enter source number"
    cat -n $date.srcs
    set srcNum=$<
    set bpTest=`head -$srcNum $date.srcs | tail -1`
    echo $bpTest
    ifrm $date.bpTst
    uvcat vis=$date.lc out=$date.bpTst select="source($bpTest)"
    gpcopy vis=$date.bp out=$date.bpTst options=nocal
    mfcal vis=$date.bpTst options=nopass flux=1,34,0
    echo ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"
    echo ">>>>> $bpTest bandpass after calibration from $bp"
    echo ">>>>>    Displayed in 30min intervals"
    echo "   Press enter to continue" ; set ans=$<
    smauvspec vis=$date.bpTst device=/xw nxy=5,5 axis=freq,both \
    interval=30 yrange=0.5,1.5
    echo "   Press enter to continue or ^C to exit without applying BP cal" 
    set ans=$<
    \rm -r $date.bpTst
endif

##################################################
BPCOPY:
echo "###>BPCOPY" >> $log

### Apply BP
rmcp $date.lc $date.bpc0
gpcopy vis=$date.bp out=$date.bpc0 options=nocal

### Regenerate wideband
ifrm $date.bpc 
uvcal vis=$date.bpc0 out=$date.bpc options=avechan
\rm -r $date.bpc0

##################################################
FLUX:
echo "###>FLUX" >> $log
Tout "### Still need a pipeline from Rudy model to this script" $log
Tout "### Rudy model variation is several percent in a day" $log
Tout "##### For now, do a very rough thing that ignores resolution" $log

## Create mars and phase cal file
ifrm $date.mars
uvaver vis=$date.bpc out=$date.mars select='source(mars,0750+125)' interval=5

####################rmcp $date.mars mars.test
####################set marsFlux=`bootflux vis=mars.test line=ch,1,1,45 taver=1000 primary=mars |grep "Calib Flux" |head -1 | awk '{print $3}'`
####################\rm  -r mars.test
set marsFlux=5.56
Tout "##### Found mars flux of $marsFlux Jy" $log
Tout "########## HARD CODED BY DPM FOR LAPTOP" $log

## gain calibrate mars data
selfcal vis=$date.mars select='source(0750+125)' axis=ph
ifrm $date.mars2
uvcal vis=$date.mars out=$date.mars2
mfcal vis=$date.mars2 options=nopass interval=2000 flux=$marsFlux,33.938,2
if !($exp == 'y') then
    gpplt vis=$date.mars2 device=/xw yrange=0.8,1.2 nxy=5,3
    echo ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"
    echo ">>>>> Plotting gains derived from mars"
    echo ">>>>>    Assumed flux of $marsFlux Jy"
    echo "   Press enter to apply to whole track" ; set ans=$<
endif

### Apply to data
rmcp $date.bpc $date.cal0
gpcopy vis=$date.mars2 out=$date.cal0 
ifrm $date.cal
uvaver vis=$date.cal0 out=$date.cal line=ch,3,1,15
\rm -r $date.cal0
echo "##### Antenna gains derived from mars:" >> $log
gplist vis=$date.mars2 >> $log

######################################################################

if !($exp == 'y') then
    echo ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"
    echo ">>>>> Plot calibrated data"
    smauvplt vis=$date.cal device=/xw line=ch,1,1,3 options=nob axis=uvd,amp \
    select='-source(jupiter)'
endif
smauvplt vis=$date.cal device=$date"_trk.ps/cps" line=ch,1,1,3 options=nob \
  select='-source(jupiter,mars)'

######################################################################
SPLIT:
echo "###>SPLIT" >> $log
if !($exp == 'y') then
    echo ">>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>>"
    echo -n ">>>>> Ready to split into single-source files?"; set ans=$<
    Tout "##### Splitting $date.cal to directory: $date" $log
endif
if !(-e $date) mkdir $date
cd $date
foreach src (`cat ../$date.srcs`)
    ifrm $src.mir
    uvcat vis=../$date.cal out=$src.mir select="source($src)"
end
cd ..
exit 0