import os, sys
from AIPS import AIPS
from AIPSTask import AIPSTask as task
from AIPSData import AIPSUVData as UV
from AIPSData import AIPSImage as IM
import time
import numpy as np

vlba_antennas = [None, 1,4,5,6,7,8,9,10,12]
#Ant   1 = BR       BX=  -219871.9213 BY= -1511722.6159 BZ=  1040436.6535
#Ant   2 = EB       BX= -5314846.1171 BY=  3909576.3105 BZ=  1214053.8483
#Ant   3 = GB       BX=   520349.7471 BY=  1635855.3622 BZ=   257352.3463
#Ant   4 = HN       BX=   -38264.6608 BY=  2118834.4173 BZ=   635929.0934
#Ant   5 = KP       BX=  1077912.6895 BY= -1190130.4420 BZ=  -329049.0013
#Ant   6 = LA       BX=   931959.8658 BY=  -660432.7647 BZ=    22746.8098
#Ant   7 = NL       BX=   516970.9535 BY=   609445.1872 BZ=   540473.9514
#Ant   8 = OV       BX=   590069.3920 BY= -1685269.9115 BZ=   152240.2949
#Ant   9 = PT       BX=  1000658.2426 BY=  -843189.7648 BZ=  -110965.2517
#Ant  10 = SC       BX=   809106.5758 BY=  3427594.6098 BZ= -1753637.4099
#Ant  11 = Y        BX=  1021321.6946 BY=  -799689.3805 BZ=  -131501.4307
#Ant  12 = MK       BX=  -895366.4616 BY= -5010871.0077 BZ= -1538079.9035


whattodo = {'load_data': True, # ->CL1
            'msort': True,
            'indxr': True,
            'preflagging': True,
            'ion_corr': True, # ->CL2
            'eops_corr': True, # ->CL3
            'par_ang_cor': True, # ->CL4
            'accor' : True, # SN1, CL5
            'mergecal': True, # HAS TO BE RUN MANUALLY
            'apcal' : True, # SN2, CL6
            # All ants have ampcal, some spikes.
            'fring_ampcal' : True, #->SN3, CL7
            'bpass_ampcal': True, # BP 1
            'image_ampcal_1': True, 
            'selfcal_ampcal' : True, # -> SN4, CL8
            'image_ampcal_2': True,
            'fring_phasecal' : True, #->SN5, CL9
            'image_phasecal_1': True,
            'selfcal_phasecal' : True, # -> SN6,7, CL10
            'image_phasecal_2': True, # 
            'image_ARP220': True,
            'selfcal_arp220': True, # ->SN8, CL11
            'image_ARP220_final': True,
            'fring_frcal' : True, #->SN9, CL12
            'image_frcal_1': True, 
            'selfcal_frcal' : True, # -> SN10, CL13
            'image_frcal_2': True, 
            'export_results': True,
            'save_calsols': True,
            }
 
########## Initialize observation data ##########
AIPS.userno = 1000
clint = 10.0/60 # Seconds
NAME = 'GC028B_U'
CLASS = 'UVDATA'
DISK = 1
AMPCALIM = NAME + '_AC'
PCALIM = NAME + '_PC'
TCALIM = NAME + '_TC'
FRCALIM = NAME + '_FC'
OUTPREFIX = NAME
target = 'ARP220'
ampcal = 'J1516+1932'
phasecal = 'J1532+2344'
frcal = 'J1613+3412'
sorted_class = 'MSORT'
logfile = 'log_'+OUTPREFIX+'_PTred.log'
refant = 1
uvrange = [None, 10000, 0]
#########^^^ END OF CONFIG VARIABLES ^^^################
tic = time.time()
########## Initialize Log file ##########
try:
 os.system('cp ' + logfile + ' ' + logfile + '.old')
 os.system('rm -f ' + logfile)
except:
 pass
AIPS.log = open(logfile, 'a')
AIPS.log.write('whattodo = '+repr(whattodo)+'\n')

########## LOAD DATA ##########
if whattodo['load_data']:
    fitld = task('fitld')
    fitld.default()
    outdata = UV(NAME, CLASS, DISK, 1)
    fitld.datain = '/data1/eskil/RAWDATA/GC028B/GC028B_UBAND.uvfits'
    fitld.outdata = outdata
    fitld.ncount = 1
    fitld.digicor = 1
    fitld.clint = clint
    fitld.wtthresh = 0.7
    fitld.doconcat = 1
    fitld.douvcomp = -1
    if outdata.exists():
        print 'WARNING: Removing existing file ' + NAME + '.' + CLASS
        outdata.zap()
        print '         File removed. Proceeding with FITLD...'
    fitld.go()

if whattodo['msort']:
    msort = task('msort')
    data = UV(NAME, CLASS, DISK, 1)
    outdata = UV(NAME, sorted_class, DISK, 1)
    msort.default()
    msort.indata = data
    msort.outdata = outdata
    msort.sort = 'TB'
    if outdata.exists():
        print 'WARNING: Removing existing file ' + NAME + '.' + sorted_class
        outdata.zap()
        print '         File removed. Proceeding with FITLD...'
    msort.go()
    if data.exists():
        print 'WARNING: Removing existing file ' + NAME + '.' + CLASS
        data.zap()
        print '         File removed. Proceeding with FITLD...'

if whattodo['indxr']:
    indxr = task('indxr')
    indxr.default()
    data = UV(NAME, sorted_class, DISK, 1)
    indxr.indata = data
    indxr.cparm[3] = clint
    indxr.go()

########## PREFLAGGING ##########
#Flag data known to be bad (whole antennas etc.)
if whattodo['preflagging']:
    data = UV(NAME, sorted_class, DISK, 1)
    infgver = 1
    outfgver = 2
    # Remove all FG tables higher than infgver.
    for i in range(data.table_highver('FG'), infgver, -1):
     data.zap_table('FG', i)
    data = UV(NAME, sorted_class, DISK, 1)
    # Copy flag table 1 to flagtable 2
    tacop = task('tacop')
    tacop.default()
    tacop.indata = data
    tacop.inext = 'FG'
    tacop.inver = infgver
    tacop.ncount = 1
    tacop.outdata = data
    tacop.outvers = outfgver
    tacop.go()
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antennas = [None, 11]
    uvflg.timeran = [None, 0, 14, 15, 0, 0, 14, 18, 0]
    uvflg.outfgver = outfgver
    uvflg.reason = 'badamp'
    uvflg.opcode = 'FLAG'
    #uvflg.go()
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antennas = [None, 2]
    uvflg.timeran = [None, 0, 14, 15, 0, 0, 14, 24, 0]
    uvflg.outfgver = outfgver
    uvflg.reason = 'badamp'
    uvflg.opcode = 'FLAG'
    #uvflg.go()
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antennas = [None, 2]
    uvflg.timeran = [None, 0, 12, 40, 0, 0, 12, 50, 0]
    uvflg.outfgver = outfgver
    uvflg.reason = 'ampspike'
    uvflg.opcode = 'FLAG'
    uvflg.go()
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antennas = [None, 3]
    uvflg.timeran = [None, 0, 20, 40, 0, 0, 20, 50, 0]
    uvflg.outfgver = outfgver
    uvflg.reason = 'ampspike'
    uvflg.opcode = 'FLAG'
    uvflg.go()
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antennas = [None, 1, 3, 5]
    uvflg.timeran = [None, 0, 9, 40, 0, 0, 10, 10, 0]
    uvflg.outfgver = outfgver
    uvflg.reason = 'ampspike'
    uvflg.opcode = 'FLAG'
    uvflg.go()
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antennas = [None, 1, 3]
    uvflg.timeran = [None, 0, 16, 50, 0, 0, 17, 10, 0]
    uvflg.outfgver = outfgver
    uvflg.reason = 'ampspike'
    uvflg.opcode = 'FLAG'
    uvflg.go()
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antennas = [None, 11]
    uvflg.timeran = [None, 0, 9, 40, 0, 0, 10, 10, 0]
    uvflg.outfgver = outfgver
    uvflg.reason = 'ampspike'
    uvflg.opcode = 'FLAG'
    uvflg.go()
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antennas = [None, 5]
    uvflg.timeran = [None, 0, 16, 00, 0, 0, 20, 0, 0]
    uvflg.outfgver = outfgver
    uvflg.reason = 'ampspike'
    uvflg.opcode = 'FLAG'
    uvflg.go()
    
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antennas = [None, 8]
    uvflg.timeran = [None, 0, 0, 0, 0, 0, 11, 0, 0]
    uvflg.outfgver = outfgver
    uvflg.reason = 'ampspike'
    uvflg.opcode = 'FLAG'
    uvflg.go()
    
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antennas = [None, 9]
    uvflg.timeran = [None, 0, 0, 0, 0, 0, 9, 30,0]
    uvflg.outfgver = outfgver
    uvflg.reason = 'ampspike'
    uvflg.opcode = 'FLAG'
    uvflg.go()
    
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antennas = [None, 11]
    uvflg.timeran = [None, 0, 9, 0, 0, 0, 10, 20, 0]
    uvflg.bif = 4
    uvflg.eif = 4
    uvflg.outfgver = outfgver
    uvflg.reason = 'ampspike'
    uvflg.opcode = 'FLAG'
    uvflg.go()
    
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antennas = [None, 11]
    uvflg.timeran = [None, 0, 20, 20, 0, 0, 22, 0, 0]
    uvflg.bif = 4
    uvflg.eif = 4
    uvflg.outfgver = outfgver
    uvflg.reason = 'ampspike'
    uvflg.opcode = 'FLAG'
    uvflg.go()
    
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antennas = [None, 2,3,11]
    uvflg.timeran = [None, 0, 14, 30, 0, 0, 14, 31, 30]
    uvflg.outfgver = outfgver
    uvflg.reason = 'ampspike'
    uvflg.opcode = 'FLAG'
    uvflg.go()
    
    # Quack
    quack = task('quack')
    quack.default()
    quack.indata = data
    quack.antenna = [None, 2,3]
    quack.opcode = 'BEG'
    quack.aparm[2] = 30.0/60
    quack.go()
    ## Quack
    quack = task('quack')
    quack.default()
    quack.indata = data
    quack.antenna = [None, 11]
    quack.opcode = 'BEG'
    quack.timeran = [None, 0, 14, 0, 0, 0, 17, 0, 0]
    quack.aparm[2] = 30.0/60
    quack.go()
    quack = task('quack')
    quack.default()
    quack.indata = data
    quack.antenna = [None, 3]
    quack.opcode = 'ENDB'
    quack.aparm[2] = 10.0/60
    quack.go()
    
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antenna = [None, 1]
    uvflg.timeran = [None, 0, 0, 0, 0, 0, 10, 0, 0]
    uvflg.outfgver = outfgver
    uvflg.reason = 'ampspike'
    uvflg.opcode = 'FLAG'
    uvflg.go()
    #uvflg = task('uvflg')
    #uvflg.default()
    #uvflg.indata = data
    #uvflg.antenna = [None, 4]
    #uvflg.timeran = [None, 0, 20, 15, 0, 0, 22, 0, 0]
    #uvflg.outfgver = outfgver
    #uvflg.reason = 'ampspike'
    #uvflg.opcode = 'FLAG'
    #uvflg.go()
    #uvflg = task('uvflg')
    #uvflg.default()
    #uvflg.indata = data
    #uvflg.antenna = [None, 7,12]
    #uvflg.timeran = [None, 0, 22, 0, 0, 0, 23, 0, 0]
    #uvflg.outfgver = outfgver
    #uvflg.reason = 'ampspike'
    #uvflg.opcode = 'FLAG'
    #uvflg.go()
    #
    ## Flag bad times antenna 6
    #uvflg = task('uvflg')
    #uvflg.default()
    #uvflg.indata = data
    #uvflg.antenna = [None, 6,10]
    #uvflg.outfgver = outfgver
    #uvflg.reason = 'ampspike'
    #uvflg.opcode = 'FLAG'
    #uvflg.timeran = [None, 0, 19, 0, 0, 0, 23, 0, 0]
    #uvflg.go()
    #uvflg = task('uvflg')
    #uvflg.default()
    #uvflg.indata = data
    #uvflg.antenna = [None, 6]
    #uvflg.outfgver = outfgver
    #uvflg.reason = 'ampspike'
    #uvflg.opcode = 'FLAG'
    #uvflg.timeran = [None, 0, 9, 45, 0, 0, 10, 0, 0]
    #uvflg.go()
    #uvflg = task('uvflg')
    #uvflg.default()
    #uvflg.indata = data
    #uvflg.antenna = [None, 10]
    #uvflg.outfgver = outfgver
    #uvflg.reason = 'ampspike'
    #uvflg.opcode = 'FLAG'
    #uvflg.timeran = [None, 0, 18, 25, 0, 0, 18, 32, 0]
    #uvflg.go()
    
    
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antenna = [None, 1]
    uvflg.outfgver = outfgver
    uvflg.reason = 'ampspike'
    uvflg.opcode = 'FLAG'
    uvflg.timeran = [None, 0, 21, 50, 0, 0, 23, 0, 0]
    #uvflg.go()
    
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antenna = [None, 6]
    uvflg.outfgver = outfgver
    uvflg.reason = 'ampspike'
    uvflg.opcode = 'FLAG'
    uvflg.timeran = [None, 0, 19, 40, 0, 0, 20, 0, 0]
    #uvflg.go()
    uvflg = task('uvflg')
    uvflg.default()
    uvflg.indata = data
    uvflg.antenna = [None, 6]
    uvflg.outfgver = outfgver
    uvflg.reason = 'ampspike'
    uvflg.opcode = 'FLAG'
    uvflg.timeran = [None, 0, 20, 40, 0, 0, 21, 0, 0]
    #uvflg.go()
    

########## IONOSPHERIC CORRECTIONS ##########
# => CL2
if whattodo['ion_corr']:
    clinver = 1
    cloutver = 2
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all CL tables higher than cloutver-1.
    for i in range(data.table_highver('CL'), cloutver-1, -1):
        data.zap_table('CL', i)
    data = UV(NAME, sorted_class, DISK, 1)
    tecor = task('tecor')
    tecor.default()
    tecor.indata = data
    tecor.gainver = clinver
    tecor.gainuse = cloutver
    # Observation was 28-DEC-2006. According to http://mistupid.com/calendar/dayofyear.htm
    # This means day of the year DDD=362
    # Fetch and uncompress file for this day
    doy = str(362)
    year = str(2006)
    tecfile = 'jplg' + doy + '0.'+year[-2:]+'i'
    if not os.path.exists(tecfile):
        os.system('wget ftp://cddis.gsfc.nasa.gov/gps/products/ionex/'+year+'/'+doy+'/' + tecfile + '.Z')
        os.system('gunzip ' + tecfile+ '.Z')
    tecor.infile = 'PWD:' + tecfile
    tecor.nfiles = 1
    tecor.go()
########## EARTH ORIENTATION PARAMETERS CORRECTION ##########
# => CL3
if whattodo['eops_corr']:
    clinver = 2
    cloutver = 3
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all CL tables higher than cloutver-1.
    for i in range(data.table_highver('CL'), cloutver-1, -1):
        data.zap_table('CL', i)
    data = UV(NAME, sorted_class, DISK, 1)
    clcor = task('clcor')
    clcor.default()
    clcor.indata = data
    clcor.gainver = clinver
    clcor.gainuse = cloutver
    clcor.opcode = 'EOPS'
    clcor.clcorprm = [None, 1, 0] 
    if not os.path.exists('./usno_finals.erp'):
        os.system('wget http://gemini.gsfc.nasa.gov/solve_save/usno_finals.erp')
    clcor.infile = 'PWD:usno_finals.erp'
    clcor.go()
 
########## PARALLACTIC ANGLE CORRECTION ##########
# => CL4
if whattodo['par_ang_cor']:
    clinver = 3
    cloutver = 4
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all CL tables higher than cloutver-1.
    for i in range(data.table_highver('CL'), cloutver-1, -1):
        data.zap_table('CL', i)
    data = UV(NAME, sorted_class, DISK, 1)
    clcor = task('clcor')
    clcor.default()
    clcor.indata = data
    clcor.gainver = clinver
    clcor.gainuse = cloutver
    clcor.opcode = 'PANG'
    clcor.clcorprm = [None, 1, 0] 
    if not os.path.exists('./usno_finals.erp'):
     os.system('wget http://gemini.gsfc.nasa.gov/solve_save/usno_finals.erp')
    clcor.infile = 'PWD:usno_finals.erp'
    clcor.go()

############# AMPLITUDE calibration ##############
# => CL5 
# => SN1(accor),
if whattodo['accor']:
    # Define initial versionnumbers, increased further down
    clinver = 4
    cloutver = 5
    snoutver = 1
## RUN ACCOR
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all SN tables higher than version snoutver-1.
    for i in range(data.table_highver('SN'), snoutver-1, -1):
        data.zap_table('SN', i)
    data = UV(NAME, sorted_class, DISK, 1)
    accor = task('accor')
    accor.indata = data
    accor.solint = 2
    accor.go()
    # Now we have a new SN-table! 
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all CL tables higher than version cloutver-1.
    for i in range(data.table_highver('CL'), cloutver-1, -1):
        data.zap_table('CL', i)
    data = UV(NAME, sorted_class, DISK, 1)
    # Make CL table
    clcal = task('clcal') 
    clcal.indata = data
    clcal.snver = snoutver
    clcal.invers = snoutver
    clcal.gainver = clinver
    clcal.gainuse = cloutver
    clcal.refant = refant
    clcal.interpol = '2PT'
    clcal.go()

## MERGECAL IS REQUIRED TO OBTAIN GC AND TY versions 2 from versions 1, else APCAL may fail.

if whattodo['mergecal']:
    ans = raw_input('Please run the procedure MERGECAL in AIPS before proceeding with apcal. Type Yes if you have run mergecal already.')
    if ans=='Yes':
        print 'Great!'
    else:
        print 'Mergecal needed! Exiting...'
        sys.exit()

if whattodo['apcal']:
## RUN APCAL 
    #redefine versions
    clinver = 5
    cloutver = 6
    snoutver = 2
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all SN tables higher than version snoutver-1.
    for i in range(data.table_highver('SN'), snoutver-1, -1):
        data.zap_table('SN', i)
    apcal = task('apcal')
    apcal.indata = data
    apcal.tyver = 2
    apcal.gcver = 2
    apcal.snver = snoutver
    apcal.invers = 1 # WX table version
    apcal.solint = 2 # Same as accor
    apcal.go()
    # Now we have a new SN-table! 
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all CL tables higher than version cloutver-1.
    for i in range(data.table_highver('CL'), cloutver-1, -1):
        data.zap_table('CL', i)
    # Make CL table
    clcal = task('clcal') 
    clcal.indata = data
    clcal.snver = snoutver
    clcal.invers = snoutver
    clcal.gainver = clinver
    clcal.gainuse = cloutver
    clcal.refant = refant
    clcal.interpol = '2PT'
    clcal.go()

#   FRING on AMPCAL
# => SN3
if whattodo['fring_ampcal']:
    clinver = 6
    cloutver = 7
    snoutver = 3
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all SN tables higher than version snoutver-1.
    for i in range(data.table_highver('SN'), snoutver-1, -1):
        data.zap_table('SN', i)
    data = UV(NAME, sorted_class, DISK, 1)
    fring=task('fring')
    fring.default()
    fring.indata = data
    fring.docalib = 1
    fring.gainuse = clinver
    fring.calsour = [None, ampcal]
    fring.refant = refant
    fring.solint = 5 # Longer than scan
    fring.snver = snoutver
    fring.go()

    # Apply SN-table
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all CL tables higher than version cloutver-1.
    for i in range(data.table_highver('CL'), cloutver-1, -1):
        data.zap_table('CL', i)
    data = UV(NAME, sorted_class, DISK, 1)
    # Make CL table
    clcal = task('clcal') 
    clcal.default()
    clcal.indata = data
    clcal.snver = snoutver
    clcal.invers = snoutver
    clcal.calsour = [None, ampcal]
    clcal.sour = [None, ''] # All sources
    clcal.gainver = clinver
    clcal.gainuse = cloutver
    clcal.refant = refant
    clcal.go()

if whattodo['bpass_ampcal']:
    clinver = 7
    bpoutver = 1
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all BP tables higher than version bpoutver-1.
    for i in range(data.table_highver('BP'), bpoutver-1, -1):
        data.zap_table('BP', i)
    data = UV(NAME, sorted_class, DISK, 1)
    data = UV(NAME, sorted_class, DISK, 1)
    bpass = task('bpass')
    bpass.default()
    bpass.indata = data
    bpass.calsour = [None,ampcal]
    bpass.timeran = [None, 0] # Use all scans to get all ants
    bpass.gainuse = clinver
    bpass.docal = 1
    bpass.refant = 1
    bpass.solint = -1 # Whole timerange
    bpass.flagver = 0
    bpass.outvers = bpoutver
    bpass.doband = -1
    bpass.bpassprm[10]=1 # Normalize amps using all channels
    bpass.go()

if whattodo['image_ampcal_1']:
    clinver = 7
    imagename = AMPCALIM + '1'
    imagedisk = 1
    imageseq = 1
    data = UV(NAME, sorted_class, DISK, 1)
    imagr = task('imagr')
    imagr.default()
    imagr.indata = data
    imagr.sources = [None, ampcal]
    imagr.gainuse = clinver
    imagr.flagver = 0
    imagr.nchav = 512
    imagr.stokes = 'I'
    imagr.cellsize = [None, 0.05e-3, 0.05e-3]
    imagr.imsize = [None, 1024,1024]
    imagr.nbox = 1
    imagr.clbox = [None, [None, 455,441,579,595]]
    imagr.docalib = 1
    imagr.dotv = -1
    imagr.outname = imagename
    imagr.outdisk = imagedisk
    imagr.outseq = imageseq
    imagr.niter = 100
    imagr.doband = 1
    imagr.bpver = 1
    imagr.uvrange = uvrange
    imagr.antenna = vlba_antennas
    imagr.baselin = imagr.antenna
    imagr.robust = 0.5
    if imagr.niter == 0:
        imageclass = 'IIM001'
    else:
        imageclass = 'ICL001'
    # Remove previous old image
    image = IM(imagename, imageclass, imagedisk, imageseq)
    if image.exists():
     image.zap()
    # Also remove old beam image
    image = IM(imagename, 'IBM001', imagedisk, imageseq)
    if image.exists():
     image.zap()
    imagr.go()


if whattodo['selfcal_ampcal']:
    clinver = 7
    cloutver = 8
    snoutver = 4
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all SN tables higher than version snoutver-1.
    for i in range(data.table_highver('SN'), snoutver-1, -1):
        data.zap_table('SN', i)
    # FIT ANTENNA 3, with proper timeran not to flag things
    data = UV(NAME, sorted_class, DISK, 1)
    calib = task('calib')
    calib.default()
    calib.indata = data
    calib.in2data = IM(AMPCALIM + '1', 'ICL001', 1, 1)
    calib.calsour = [None, ampcal]
    calib.refant = refant
    calib.docalib = 1
    calib.gainuse = clinver
    calib.flagver = 0
    calib.snver = snoutver
    calib.solmode = 'A&P'
    calib.solint = 5 # 
    calib.aparm[7] = 5 # SNR
    calib.cmethod = 'DFT'
    calib.aparm[6] = 1 # Print level
    calib.dofit = [None, 2,3,11] # Fit 
    calib.doband = 1
    calib.bpver = 1
    calib.go()
    # Apply SN-table
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all CL tables higher than version cloutver-1.
    for i in range(data.table_highver('CL'), cloutver-1, -1):
        data.zap_table('CL', i)
    data = UV(NAME, sorted_class, DISK, 1)
    # Make CL table
    clcal = task('clcal') 
    clcal.default()
    clcal.indata = data
    clcal.snver = snoutver
    clcal.invers = snoutver
    clcal.antenna = [None, 2,3,11]
    clcal.calsour = [None, ampcal] # Use sols from Ampcal 
    clcal.sour = [None, ''] # Apply to all sources
    clcal.gainver = clinver
    clcal.gainuse = cloutver
    clcal.refant = refant
    clcal.go()

if whattodo['image_ampcal_2']:
    clinver = 8
    imagename = AMPCALIM+'2'
    imagedisk = 1
    imageseq = 1
    data = UV(NAME, sorted_class, DISK, 1)
    imagr = task('imagr')
    imagr.default()
    imagr.indata = data
    imagr.sources = [None, ampcal]
    imagr.gainuse = clinver
    imagr.flagver = 0
    imagr.nchav = 512
    imagr.stokes = 'I'
    imagr.cellsize = [None, 0.05e-3, 0.05e-3]
    imagr.imsize = [None, 1024,1024]
    imagr.docalib = 1
    imagr.dotv = -1
    imagr.doband = 1
    imagr.bpver = 1
    imagr.outname = imagename
    imagr.outdisk = imagedisk
    imagr.outseq = imageseq
    imagr.niter = 300
    imagr.uvrange = uvrange
    imagr.robust = 0.5
    if imagr.niter == 0:
        imageclass = 'IIM001'
    else:
        imageclass = 'ICL001'
    # Remove previous old image
    image = IM(imagename, imageclass, imagedisk, imageseq)
    if image.exists():
     image.zap()
    # Also remove old beam image
    image = IM(imagename, 'IBM001', imagedisk, imageseq)
    if image.exists():
     image.zap()
    imagr.go()

#   FRING on PHASECAL
if whattodo['fring_phasecal']:
    clinver = 8 # With ampcal corrections
    cloutver = 9
    snoutver = 5
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all SN tables higher than version snoutver-1.
    for i in range(data.table_highver('SN'), snoutver-1, -1):
        data.zap_table('SN', i)
    data = UV(NAME, sorted_class, DISK, 1)
    fring=task('fring')
    fring.default()
    fring.indata = data
    fring.docalib = 1
    fring.gainuse = clinver
    fring.calsour = [None, phasecal]
    fring.refant = refant
    fring.solint = 3 # Longer than scan
    fring.dparm[2] = 150#ns, delay win, already fitted bulk on ampcal
    fring.dparm[3] = 50#mHz, rate win
    #fring.aparm[3] = 1 # AVG RR,LL
    #fring.aparm[4] = 1 # AVG freqs in IFs
    fring.aparm[5] = 1 # AVG IFs
    fring.snver = snoutver
    fring.doband = 1
    fring.bpver = 1
    fring.go()

    # Apply SN-table
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all CL tables higher than version cloutver-1.
    for i in range(data.table_highver('CL'), cloutver-1, -1):
        print "WARNING: Deleting CL" + str(i)
        data.zap_table('CL', i)
    data = UV(NAME, sorted_class, DISK, 1)
    # Make CL table
    clcal = task('clcal') 
    clcal.default()
    clcal.indata = data
    clcal.snver = snoutver
    clcal.invers = snoutver
    clcal.calsour = [None, phasecal]
    clcal.sour = [None, phasecal, target]
    clcal.gainver = clinver
    clcal.gainuse = cloutver
    clcal.refant = refant
    clcal.go()

if whattodo['image_phasecal_1']:
    clinver = 9
    imagename = PCALIM + '1'
    imagedisk = 1
    imageseq = 1
    data = UV(NAME, sorted_class, DISK, 1)
    imagr = task('imagr')
    imagr.default()
    imagr.indata = data
    imagr.sources = [None, phasecal]
    imagr.gainuse = clinver
    imagr.flagver = 0
    imagr.nchav = 512
    imagr.stokes = 'I'
    imagr.cellsize = [None, 0.05e-3, 0.05e-3]
    imagr.imsize = [None, 1024,1024]
    imagr.docalib = 1
    imagr.dotv = -1
    imagr.doband = 1
    imagr.bpver = 1
    imagr.outname = imagename
    imagr.outdisk = imagedisk
    imagr.outseq = imageseq
    imagr.niter = 300
    imagr.uvrang = uvrange
    imagr.robust = 0.5
    if imagr.niter == 0:
        imageclass = 'IIM001'
    else:
        imageclass = 'ICL001'
    # Remove previous old image
    image = IM(imagename, imageclass, imagedisk, imageseq)
    if image.exists():
     image.zap()
    # Also remove old beam image
    image = IM(imagename, 'IBM001', imagedisk, imageseq)
    if image.exists():
     image.zap()
    imagr.go()


if whattodo['selfcal_phasecal']:
    clinver = 9
    cloutver = 10
    snoutver = 6 # ,7
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all SN tables higher than version snoutver-1.
    for i in range(data.table_highver('SN'), snoutver-1, -1):
        print "WARNING: Deleting SN" + str(i)
        data.zap_table('SN', i)
    data = UV(NAME, sorted_class, DISK, 1)
    calib = task('calib')
    calib.default()
    calib.indata = data
    calib.in2data = IM(PCALIM + '1', 'ICL001', 1, 1)
    calib.calsour = [None, phasecal]
    calib.refant = refant
    calib.docalib = 1
    calib.gainuse = clinver
    calib.flagver = 0
    calib.snver = snoutver
    calib.solmode = 'A&P'
    calib.solint = 3 # Longer than scan
    calib.aparm[7] = 4 # SNR
    calib.aparm[3] = 1 # AVG RR,LL, 
    calib.aparm[6] = 1 # Print level
    calib.aparm[5] = 1 # AVG ifs
    calib.weightit = 1
    calib.soltype = 'L1R'
    calib.doband = 1
    calib.bpver = 1
    calib.go()
    
    snsmo = task('snsmo')
    snsmo.indata = data
    snsmo.inver = snoutver
    snsmo.outver = snoutver+1
    snsmo.smotype = 'AMPL'
    snsmo.cparm[1] = 60.0/60 #in hours
    snsmo.cparm[6] = 0.2 # Max amp deviation
    snsmo.bparm[1] = 120.0/60 #in hours
    snsmo.smotype = 'GAUS'
    snsmo.doblank = -1 # Since cannot be used with AMPL
    snsmo.go()
    
    # Apply SN-table
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all CL tables higher than version cloutver-1.
    for i in range(data.table_highver('CL'), cloutver-1, -1):
        print "WARNING: Deleting CL" + str(i)
        data.zap_table('CL', i)
    data = UV(NAME, sorted_class, DISK, 1)
    # Make CL table
    clcal = task('clcal') 
    clcal.default()
    clcal.indata = data
    clcal.snver = snoutver+1
    clcal.invers = snoutver+1
    clcal.calsour = [None, phasecal]
    clcal.sour = [None, phasecal, target] # Apply to all sources
    clcal.gainver = clinver
    clcal.gainuse = cloutver
    clcal.refant = refant
    clcal.go()

if whattodo['image_phasecal_2']:
    clinver = 10
    imagename = PCALIM + '2'
    imagedisk = 1
    imageseq = 1
    data = UV(NAME, sorted_class, DISK, 1)
    imagr = task('imagr')
    imagr.default()
    imagr.indata = data
    imagr.sources = [None, phasecal]
    imagr.gainuse = clinver
    imagr.flagver = 0
    imagr.nchav = 512
    imagr.stokes = 'I'
    imagr.cellsize = [None, 0.05e-3, 0.05e-3]
    imagr.imsize = [None, 1024,1024]
    imagr.docalib = 1
    imagr.dotv = -1
    imagr.doband = 1
    imagr.bpver = 1
    imagr.outname = imagename
    imagr.outdisk = imagedisk
    imagr.outseq = imageseq
    imagr.niter = 300
    imagr.uvran = uvrange
    imagr.robust = 0.5
    if imagr.niter == 0:
        imageclass = 'IIM001'
    else:
        imageclass = 'ICL001'
    # Remove previous old image
    image = IM(imagename, imageclass, imagedisk, imageseq)
    if image.exists():
     image.zap()
    # Also remove old beam image
    image = IM(imagename, 'IBM001', imagedisk, imageseq)
    if image.exists():
     image.zap()
    imagr.go()

############### IMAGE ARP220 ######################
# Load these outside if to be available also for export later
arpfilename = TCALIM
arpfiledisk = 1
arpfileseq = 1
arpfileclass = 'ICL001'
arpfileclass2 = 'ICL002'
if whattodo['image_ARP220']:
    # FIRST CLEAN WITH SMALL BOXES
    # Remove previous old images
    image = IM(arpfilename, arpfileclass, arpfiledisk, arpfileseq)
    if image.exists():
        image.zap()
    image = IM(arpfilename, arpfileclass2, arpfiledisk, arpfileseq)
    if image.exists():
        image.zap()
    # Also remove old beam images
    image = IM(arpfilename, 'IBM001', arpfiledisk, arpfileseq)
    if image.exists():
        image.zap()
    image = IM(arpfilename, 'IBM002', arpfiledisk, arpfileseq)
    if image.exists():
        image.zap()
    # Clean with latest gainuse
    data = UV(NAME, sorted_class, DISK, 1)
    imagr = task('imagr')
    imagr.default()
    imagr.indata = data
    imagr.sources = [None, target]
    imagr.im2parm = [None, 1, 5.0, 10.0, 0] # Require island to have SNR 5 and peak SNR10. To clean (for later removal) of the strongest sources
    imagr.gainuse = 0
    imagr.flagver = 0
    imagr.nchav = 512
    imagr.doband = 1
    imagr.bpver = 1
    imagr.imsize = [None, 8192,4096]
    imagr.cellsize = [None, 0.05e-3, 0.05e-3]
    imagr.stokes = 'I'
    imagr.nfield = 2
    imagr.overlap = 1
    imagr.do3dimag = 1
    imagr.rashift = [None, 0.6, -0.35]
    imagr.decshift = [None, 0.02, 0.18]
    imagr.docalib = 1
    imagr.dotv = -1
    imagr.outname = arpfilename
    imagr.outdisk = arpfiledisk
    imagr.outseq = arpfileseq
    imagr.niter = 500
    imagr.gain = 0.1
    imagr.uvrange = uvrange # remove large-scale disturbances
    imagr.robust = 0.5
    imagr.go()
    # CLEAN WITH LARGE (NO) BOXES
    arpboxfile = ''
    data = UV(NAME, sorted_class, DISK, 1)
    imagr = task('imagr')
    imagr.default()
    imagr.indata = data
    imagr.sources = [None, target]
    imagr.gainuse = 0
    imagr.flagver = 0
    imagr.nchav = 512
    imagr.doband = 1
    imagr.bpver = 1
    imagr.imsize = [None, 8192,4096]
    imagr.cellsize = [None, 0.05e-3, 0.05e-3]
    imagr.stokes = 'I'
    imagr.nfield = 2
    imagr.overlap = 1
    imagr.do3dimag = 1
    imagr.rashift = [None, 0.6, -0.35]
    imagr.decshift = [None, 0.02, 0.18]
    imagr.docalib = 1
    imagr.dotv = -1
    imagr.outname = arpfilename
    imagr.outdisk = arpfiledisk
    imagr.outseq = arpfileseq
    imagr.outver = 1
    imagr.niter = 1000
    imagr.gain = 0.1
    imagr.bcomp = [None, 100000, 100000]
    imagr.uvrange = uvrange # remove large-scale disturbances
    imagr.robust = 0.5
    imagr.go()

if whattodo['selfcal_arp220']:
    clinver = 10
    cloutver = 11
    snoutver = 8
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all SN tables higher than version snoutver-1.
    for i in range(data.table_highver('SN'), snoutver-1, -1):
        print "WARNING: Deleting SN" + str(i)
        data.zap_table('SN', i)
    data = UV(NAME, sorted_class, DISK, 1)
    calib = task('calib')
    calib.default()
    calib.indata = data
    # Do not use any model, just a point works fine for phasecal
    calib.in2data = IM(TCALIM, 'ICL001', 1, 1)
    calib.calsour = [None, target]
    calib.nmaps = 2
    calib.refant = refant
    calib.docalib = 1
    calib.gainuse = clinver
    calib.flagver = 0
    calib.snver = snoutver
    calib.solmode = 'P'
    calib.solint = 10 # Longer than scan
    calib.aparm[1] = 3 # min No. of antennas to include first 45 min.
    calib.aparm[7] = 2 # SNR
    calib.aparm[3] = 1 # AVG RR,LL, 
    calib.aparm[6] = 1 # Print level
    calib.aparm[5] = 1 # AVG ifs
    calib.cmethod = 'DFT'
    calib.weightit = 1
    calib.soltype = 'L1R'
    calib.doband = 1
    calib.bpver = 1
    calib.uvrange = uvrange
    calib.go()
    
    # Apply SN-table
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all CL tables higher than version cloutver-1.
    for i in range(data.table_highver('CL'), cloutver-1, -1):
        print "WARNING: Deleting CL" + str(i)
        data.zap_table('CL', i)
    data = UV(NAME, sorted_class, DISK, 1)
    # Make CL table
    clcal = task('clcal') 
    clcal.default()
    clcal.indata = data
    clcal.snver = snoutver
    clcal.invers = snoutver
    clcal.calsour = [None, target]
    clcal.sour = [None, target] # Apply to target
    clcal.gainver = clinver
    clcal.gainuse = cloutver
    clcal.refant = refant
    clcal.go()

############### IMAGE ARP220 ######################
# Load these outside if to be available also for export later
arpfilename = OUTPREFIX + '_FIN'
arpfiledisk = 1
arpfileseq = 1
arpfileclass = 'ICL001'
arpfileclass2 = 'ICL002'
if whattodo['image_ARP220_final']:
    # FIRST CLEAN WITH SMALL BOXES
    # Remove previous old images
    image = IM(arpfilename, arpfileclass, arpfiledisk, arpfileseq)
    if image.exists():
        image.zap()
    image = IM(arpfilename, arpfileclass2, arpfiledisk, arpfileseq)
    if image.exists():
        image.zap()
    # Also remove old beam images
    image = IM(arpfilename, 'IBM001', arpfiledisk, arpfileseq)
    if image.exists():
        image.zap()
    image = IM(arpfilename, 'IBM002', arpfiledisk, arpfileseq)
    if image.exists():
        image.zap()
    # Clean with latest gainuse
    data = UV(NAME, sorted_class, DISK, 1)
    imagr = task('imagr')
    imagr.default()
    imagr.indata = data
    imagr.sources = [None, target]
    imagr.im2parm = [None, 1, 5.0, 10.0, 0] # Require island to have SNR 5 and peak SNR10. To clean (for later removal) of the strongest sources
    imagr.gainuse = 0
    imagr.flagver = 0
    imagr.nchav = 512
    imagr.doband = 1
    imagr.bpver = 1
    imagr.imsize = [None, 8192, 8192]
    imagr.cellsize = [None, 0.05e-3, 0.05e-3]
    imagr.stokes = 'I'
    imagr.nfield = 2
    imagr.overlap = 1
    imagr.do3dimag = 1
    imagr.rashift = [None, 0.6, -0.35]
    imagr.decshift = [None, 0.02, 0.18]
    imagr.docalib = 1
    imagr.dotv = -1
    imagr.outname = arpfilename
    imagr.outdisk = arpfiledisk
    imagr.outseq = arpfileseq
    imagr.niter = 500
    imagr.gain = 0.1
    imagr.uvrange = uvrange  # remove large-scale disturbances
    imagr.robust = 5
    imagr.go()
    # CLEAN WITH LARGE (NO) BOXES
    arpboxfile = ''
    data = UV(NAME, sorted_class, DISK, 1)
    imagr = task('imagr')
    imagr.default()
    imagr.indata = data
    imagr.sources = [None, target]
    imagr.gainuse = 0
    imagr.flagver = 0
    imagr.nchav = 512
    imagr.doband = 1
    imagr.bpver = 1
    imagr.imsize = [None, 8192, 8192]
    imagr.cellsize = [None, 0.05e-3, 0.05e-3]
    imagr.stokes = 'I'
    imagr.nfield = 2
    imagr.overlap = 1
    imagr.do3dimag = 1
    imagr.rashift = [None, 0.6, -0.35]
    imagr.decshift = [None, 0.02, 0.18]
    imagr.docalib = 1
    imagr.dotv = -1
    imagr.outname = arpfilename
    imagr.outdisk = arpfiledisk
    imagr.outseq = arpfileseq
    imagr.outver = 1
    imagr.niter = 1000
    imagr.gain = 0.1
    imagr.bcomp = [None, 100000, 100000]
    imagr.uvrange = uvrange # remove large-scale disturbances
    imagr.robust = 5
    imagr.go()

#   FRING on FRCAL
if whattodo['fring_frcal']:
    clinver = 11
    cloutver = 12
    snoutver = 9
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all SN tables higher than version snoutver-1.
    for i in range(data.table_highver('SN'), snoutver-1, -1):
        data.zap_table('SN', i)
    data = UV(NAME, sorted_class, DISK, 1)
    fring=task('fring')
    fring.default()
    fring.indata = data
    fring.docalib = 1
    fring.gainuse = clinver
    fring.calsour = [None, frcal]
    fring.refant = refant
    fring.solint = 10 # Longer than scan
    fring.snver = snoutver
    fring.go()

    # Apply SN-table
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all CL tables higher than version cloutver-1.
    for i in range(data.table_highver('CL'), cloutver-1, -1):
        data.zap_table('CL', i)
    data = UV(NAME, sorted_class, DISK, 1)
    # Make CL table
    clcal = task('clcal') 
    clcal.default()
    clcal.indata = data
    clcal.snver = snoutver
    clcal.invers = snoutver
    clcal.calsour = [None, frcal]
    clcal.sour = [None, frcal]
    clcal.gainver = clinver
    clcal.gainuse = cloutver
    clcal.refant = refant
    clcal.go()

if whattodo['image_frcal_1']:
    clinver = 12
    imagename = FRCALIM + '1'
    imagedisk = 1
    imageseq = 1
    data = UV(NAME, sorted_class, DISK, 1)
    imagr = task('imagr')
    imagr.default()
    imagr.indata = data
    imagr.sources = [None, frcal]
    imagr.gainuse = clinver
    imagr.flagver = 0
    imagr.nchav = 512
    imagr.stokes = 'I'
    imagr.cellsize = [None, 0.05e-3, 0.05e-3]
    imagr.imsize = [None, 1024,1024]
    imagr.docalib = 1
    imagr.dotv = -1
    imagr.outname = imagename
    imagr.outdisk = imagedisk
    imagr.outseq = imageseq
    imagr.niter = 200
    imagr.doband = 1
    imagr.bpver = 1
    imagr.uvrange = uvrange
    imagr.robust = 0.5
    imagr.nbox = 1
    imagr.clbox = [None, [None, 455,441,579,595]]
    if imagr.niter == 0:
        imageclass = 'IIM001'
    else:
        imageclass = 'ICL001'
    # Remove previous old image
    image = IM(imagename, imageclass, imagedisk, imageseq)
    if image.exists():
     image.zap()
    # Also remove old beam image
    image = IM(imagename, 'IBM001', imagedisk, imageseq)
    if image.exists():
     image.zap()
    imagr.go()

if whattodo['selfcal_frcal']:
    clinver = 12
    cloutver = 13
    snoutver = 10
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all SN tables higher than version snoutver-1.
    for i in range(data.table_highver('SN'), snoutver-1, -1):
        data.zap_table('SN', i)
    # FIT ANTENNA 3, with proper timeran not to flag things
    data = UV(NAME, sorted_class, DISK, 1)
    calib = task('calib')
    calib.default()
    calib.indata = data
    calib.in2data = IM(FRCALIM + '1', 'ICL001', 1, 1)
    calib.calsour = [None, frcal]
    calib.refant = refant
    calib.docalib = 1
    calib.gainuse = clinver
    calib.flagver = 0
    calib.snver = snoutver
    calib.solmode = 'A&P'
    calib.solint = 0.5 # 
    calib.aparm[7] = 5 # SNR
    calib.cmethod = 'DFT'
    calib.aparm[6] = 1 # Print level
    calib.doband = 1
    calib.bpver = 1
    calib.go()
    # Apply SN-table
    data = UV(NAME, sorted_class, DISK, 1)
    # Remove all CL tables higher than version cloutver-1.
    for i in range(data.table_highver('CL'), cloutver-1, -1):
        data.zap_table('CL', i)
    data = UV(NAME, sorted_class, DISK, 1)
    # Make CL table
    clcal = task('clcal') 
    clcal.default()
    clcal.indata = data
    clcal.snver = snoutver
    clcal.invers = snoutver
    clcal.calsour = [None, frcal]
    clcal.sour = [None, frcal]
    clcal.gainver = clinver
    clcal.gainuse = cloutver
    clcal.refant = refant
    clcal.go()

if whattodo['image_frcal_2']:
    clinver = 13
    imagename = FRCALIM+'2'
    imagedisk = 1
    imageseq = 1
    data = UV(NAME, sorted_class, DISK, 1)
    imagr = task('imagr')
    imagr.default()
    imagr.indata = data
    imagr.sources = [None, frcal]
    imagr.gainuse = clinver
    imagr.flagver = 0
    imagr.nchav = 512
    imagr.stokes = 'I'
    imagr.cellsize = [None, 0.05e-3, 0.05e-3]
    imagr.imsize = [None, 1024,1024]
    imagr.docalib = 1
    imagr.dotv = -1
    imagr.doband = 1
    imagr.bpver = 1
    imagr.outname = imagename
    imagr.outdisk = imagedisk
    imagr.outseq = imageseq
    imagr.niter = 200
    imagr.uvrange = uvrange
    imagr.robust = 0.5
    imagr.nbox = 1
    imagr.clbox = [None, [None, 455,441,579,595]]
    if imagr.niter == 0:
        imageclass = 'IIM001'
    else:
        imageclass = 'ICL001'
    # Remove previous old image
    image = IM(imagename, imageclass, imagedisk, imageseq)
    if image.exists():
     image.zap()
    # Also remove old beam image
    image = IM(imagename, 'IBM001', imagedisk, imageseq)
    if image.exists():
     image.zap()
    imagr.go()

############### EXPORT RESULTS ######################
if whattodo['export_results']:
    east = IM(arpfilename, arpfileclass, arpfiledisk, arpfileseq)
    west = IM(arpfilename, arpfileclass2, arpfiledisk, arpfileseq)
   
    # Export EAST
    fittp = task('fittp')
    fittp.default()
    fittp.indata = east
    fittp.dataout = 'PWD:'+OUTPREFIX+'_EAST_IMAGE.FITS'
    if os.path.exists(fittp.dataout[4:]):
     os.system('rm ' + fittp.dataout[4:])
    fittp.go()
    # Export WEST
    fittp = task('fittp')
    fittp.default()
    fittp.indata = west
    fittp.dataout = 'PWD:'+OUTPREFIX+'_WEST_IMAGE.FITS'
    if os.path.exists(fittp.dataout[4:]):
     os.system('rm ' + fittp.dataout[4:])
    fittp.go()
    # Export FRCAL
    im = IM(NAME+'_FC2', 'ICL001', 1, 1)
    fittp = task('fittp')
    fittp.default()
    fittp.indata = im
    fittp.dataout = 'PWD:'+NAME+'_J1613+3412_IMAGE.FITS'
    if os.path.exists(fittp.dataout[4:]):
     os.system('rm ' + fittp.dataout[4:])
    fittp.go()
    # Export ACAL
    im = IM(NAME+'_AC2', 'ICL001', 1, 1)
    fittp = task('fittp')
    fittp.default()
    fittp.indata = im
    fittp.dataout = 'PWD:'+NAME+'_J1516+1932_IMAGE.FITS'
    if os.path.exists(fittp.dataout[4:]):
     os.system('rm ' + fittp.dataout[4:])
    fittp.go()
    # Export PCAL
    im = IM(NAME+'_PC2', 'ICL001', 1, 1)
    fittp = task('fittp')
    fittp.default()
    fittp.indata = im
    fittp.dataout = 'PWD:'+NAME+'_J1532+2344_IMAGE.FITS'
    if os.path.exists(fittp.dataout[4:]):
     os.system('rm ' + fittp.dataout[4:])
    fittp.go()
    
    fittp = task('fittp')
    fittp.default()
    fittp.indata = IM(NAME+'_FIN','IBM001', 1, 1)
    fittp.dataout = 'PWD:'+NAME+'_EAST_BEAM.FITS'
    if os.path.exists(fittp.dataout[4:]):
     os.system('rm ' + fittp.dataout[4:])
    fittp.go()
    fittp = task('fittp')
    fittp.default()
    fittp.indata = IM(NAME+'_FIN','IBM002', 1, 1)
    fittp.dataout = 'PWD:'+NAME+'_WEST_BEAM.FITS'
    if os.path.exists(fittp.dataout[4:]):
     os.system('rm ' + fittp.dataout[4:])
    fittp.go()

if whattodo['save_calsols']:
    snver = 8
    # Remove all PL tables
    data = UV(NAME, sorted_class, DISK, 1)

    data.zap_table('PL', -1)
    snplt = task('snplt')
    snplt.indata = data
    snplt.bif = 1
    snplt.eif = 1
    snplt.stokes = 'RR'
    snplt.optype = 'PHAS'
    snplt.sour = [None, target]
    snplt.inext = 'SN'
    snplt.inver = snver
    snplt.nplots = 6
    #snplt.pixrange = [None, -180,180]
    snplt.dotv = -1
    snplt.go()
    outname = NAME+'_ARP220_SELFCAL_SOLS.EPS'
    if os.path.exists(outname):
        os.system('rm ' + outname)
    lwpla = task('lwpla')
    lwpla.indata = data
    lwpla.plver = 1
    lwpla.inver = 99
    lwpla.outfile = 'PWD:'+outname
    #lwpla.docolor = 1
    lwpla.go()


tac = time.time()
processtime = tac-tic
print 'OVERALL PROCESS TIME: ',processtime,' SECONDS. or ', processtime/60.0, 'minutes'