#! /bin/tcsh 
#
# JAC
# 1.0 091006 - original

set version=1.0

# nice +19

cat <<EOF 

--------------- `echo $0:t | tr '[a-z]' '[A-Z]'` version $version `date +%x` -----------------

EOF

if ( $#argv < 8) then
 cat <<EOF 

Script removes bright sources from a given event file, using a source
file as the basis for this extraction, and refills these regions by
ghosting events from the surrounding area in DETX and DETY
coordinates. An extraction radius of about 25 to 40 arcseconds is
recommended.  Greater radii may be used, but this may currently result
in under-filling. At the end of the script, the ghosted events are
recast using the SAS task attcalc. The aim is to be left with a
reasonably uniform event file

This script uses the XMM-Newton SAS tasks: region evselect attcalc
 (http://xmm.vilspa.esa.es/external/xmm_sw_cal/sas_frame.shtml)
SAS version 7.0.0 is recommended.

The source file should be of a format whereby the PN, MOS1 and MOS2
information is combined (this often contains the string OMSRLI).

The SAS_CCFPATH should be correctly set.

This script uses FTOOLS.

This script uses IDL. The procedure ghostholes.pro is required. The
argument idlcodedir lets the user chose the location of this
script. IDL must be correctly set up.
(http://www.ittvis.com/idl/index.asp)

This script makes no check if the source and event file correspond to
the same observation.

The event files and source files are likely to contain the following
strings:

events - *EVLI*
source - *OMSRLI*

Calling sequence

    ./ghostholes_ind events_file source_file cal_file radius idlcodedir finalfile keepholes

Required arguments 
    events_file - event file to be used from 2XMM catalogue
    source_file - source file to be used from 2XMM catalogue
    cal_file    - calibration index file
    radius      - radius of extraction region in arcseconds
    idlcodedir  - where the idl ghostholes.pro code can be found
    clean       - cleans event files in the chip gaps or outside detector areas, Y/N
    finalfile   - name for final file. If N then default used
    keepholes   - Y/N, keep event files with bright sources removed

Advice:
    For directory names, such as idlcodedir a '/' is needed, i.e. /home/software/idl/
    The IDL astrolib is required for this code, please see
	http://idlastro.gsfc.nasa.gov/contents.html
    
---------------------------------------------------------------

Example: ghostholes, radius 25 arcseconds, clean and keep source removed events file

e.g. ghostholes_ind myevents.fits mysource.fits mycalfile.fits 25 /home/code/idl/ Y ghosted.fits Y

---------------------------------------------------------------
EOF
 exit
endif


set evlist=$1
set source=$2
set ccf=$3
set srcrad=$4
set idldir=$5
set clean=$6
set finalfile=$7
set keep=$8

set prog=$idldir'ghostit_ind.pro'

#echo $keep
#exit

if ($finalfile == 'N' || $finalfile == 'n' || $finalfile == 'No' || $finalfile == 'no') then     
    set finalfile='filled.events.fits'
endif

if ($keep == 'Y' || $keep == 'yes' || $keep == 'y' || $keep != '') then     
    set keep='Y'
endif

if ($finalfile == '') then 
    set finalfile='filled.events.fits'
endif

#-----------------------------------------#
#  Check input 
#-----------------------------------------#

if (! -e $evlist) then 
 echo "ERROR: eventlist " $evlist " does not exist...exiting..."
 exit
endif

if (! -e $source) then 
 echo "ERROR: source file " $source " does not exist...exiting..."
 exit
endif

if (! -e $ccf) then
    echo 'ERROR: ccf file " $ccf " does not exist... exiting...'
    exit
endif

if (! -e $prog) then
    echo 'ERROR: idl dir/prog " $prog " does not exist... exiting...'
    exit
endif


#-----------------------------------------#
#  Write batch file 
#-----------------------------------------#   

set myfirstBATCH_file=tempBATCH

if ( -e $myfirstBATCH_file ) rm -f $myfirstBATCH_file
touch $myfirstBATCH_file
echo 'Program (IDL) to use: ' $prog
echo .compile $prog >> $myfirstBATCH_file
# echo ghostholes_ind >> $myfirstBATCH_file
echo simple_ghostit_ind >> $myfirstBATCH_file
echo exit >> $myfirstBATCH_file

#-----------------------------------------#
# Check validity of files
#-----------------------------------------#

if ( -e VERIFY ) rm -f VERIFY
fverify $source outfile=VERIFY clobber=yes 
set valid=`cat VERIFY | grep SRCLIST | awk '{print $1}'`
set valid=`echo $valid | cut -d " " -f1`

if ( $valid != "SRCLIST" ) then
    echo 'ERROR: Source file not valid source fits'
    exit
endif

set valid=""

if ( -e VERIFY ) rm -f VERIFY
fverify $evlist outfile=VERIFY clobber=yes 
set valid=`cat VERIFY | grep EVENTS | awk '{print $1}'`
set valid=`echo $valid | cut -d " " -f1`
set nexti=`cat VERIFY | grep Header | awk '{print $1}'`

if ( $valid != "EVENTS" ) then
    echo 'ERROR: Events file not valid events fits'
    exit
endif

set valid=""

if ( -e tempfd ) rm -f tempfd
fdump $ccf'+1' tempfd 1 1-1 clobber=yes  
set valid=`grep "CALINDEX" < tempfd | head -1 | awk '{print $3}'` 
set valid=`echo $valid| cut -d "'" -f2`

if ( $valid != "CALINDEX") then
    echo 'ERROR: Calibration index file not valid calibration fits'
    exit
endif

#-----------------------------------------#
#  Extract info from evlist
#-----------------------------------------#

setenv SAS_CCF $ccf
echo 'SAS variables are set up as follows...'
printenv | grep SAS

if ( -e tempfd ) rm -f tempfd
fdump $evlist'+1' tempfd 1 1-1 clobber=yes    

set datamodeo=`grep "DATAMODE" < tempfd | head -1 | awk '{print $2}'` 
set datamode=`echo $datamodeo | cut -d "'" -f2`

set smode=`grep "SUBMODE " < tempfd | head -1 | awk '{print $3}'`
set smode=`echo $smode | cut -d "'" -f2`

    if ( $smode == 'PrimeFullWindow') set smode='FF'
    if ( $smode == 'PrimeFullWindowExt') set smode='EF'
    if ( $smode == 'PrimeFullWindowExtended') set smode='EF'

set filtero=`grep "FILTER " < tempfd | head -1 | awk '{print $3}'`
set filter=`echo $filtero | cut -d "'" -f2`
set filterk=$filter

    if ( $filter == 'Thin1' || $filter == 'Thin2' ) set filter='T'
    if ( $filter == 'Medium') set filter='M'
    if ( $filter == 'Thick') set filter='K'

set instr=`grep "INSTRUME" < tempfd | head -1 | awk '{print $2}'`
set instr=`echo $instr | cut -d "'" -f2`
if ( $instr == 'EPN' ) set instr="PN"
if ( $instr == 'EMOS1' ) set instr="M1"
if ( $instr == 'EMOS2' ) set instr="M2"

set ra_pnt = `grep "RA_PNT" < tempfd | head -1 | awk '{print $3}'` 
set dec_pnt = `grep "DEC_PNT" < tempfd | head -1 | awk '{print $3}'` 
set pa_pnt = `grep "PA_PNT" < tempfd | head -1 | awk '{print $3}'`  
set ranom=`grep REFXCRVL < tempfd | head -1 | awk '{print $2}'`
set decnom=`grep REFYCRVL < tempfd | head -1 | awk '{print $2}'`

#-----------------------------------------#
# Extract region from source file
#-----------------------------------------#

set evlist_2=tempEvents

region -w 0 -V 0 eventset=$evlist operationstyle=global radiusstyle=userfixed fixedradius=$srcrad srclisttab=$source expression="ID_INST==0" outunit="detxy" bkgregionset='region.reg'
	    
# Region file here set to that in DETX, DETY
evselect -w 0 -V 0 table=$evlist expression="region(region.reg,DETX,DETY)" filteredset=$evlist_2 updateexposure=N writedss=N withfilteredset=Y destruct=Y keepfilteroutput=T

#-----------------------------------------#
# Run IDL code
#-----------------------------------------#
if (-e preidl.fits) rm -f preidl.fits
cp tempEvents preidl.fits

set regionFile='region.reg'

set IDLcallfile=$myfirstBATCH_file
setenv IDL_STARTUP $IDLcallfile
setenv IDL_ARG1 $instr
setenv IDL_ARG2 $evlist_2
setenv IDL_ARG3 $regionFile
setenv IDL_ARG4 'tempfill.fits'

echo 'Going to IDL code'
idl 
setenv IDL_STARTUP

#-----------------------------------------#
# Replace lost header keywords
#-----------------------------------------#

fparkey "'"$datamode"'" tempfill.fits+0 DATAMODE comm='Instrument mode (IMAGING, TIMING, BURST, etc.)' add=yes
fparkey "'"$filterk"'" tempfill.fits+0 FILTER comm='Filter ID' add=yes

#-----------------------------------------#
# Cleaning
#-----------------------------------------#

if (-e preallcleaning.fits) rm -f preallcleaning.fits
cp tempfill.fits preallcleaning.fits

echo 'Go to clean ' $clean

if ( $clean == 'Y' || $clean == 'Yes' || $clean == 'yes' || $clean == 'y' ) then
    attcalc -V 0 -w 0 eventset=tempfill.fits attitudelabel=fixed refpointlabel=user withatthkset=N fixedra=0 fixeddec=0 fixedposangle=0 nominalra=0 nominaldec=0
    if ( $instr == 'PN' ) then 
	echo "PN event cleanup"	    
	evselect -w 0 -V 0 table=tempfill.fits filteredset=tempfill2.fits updateexposure=Y writedss=N withfilteredset=Y destruct=Y keepfilteroutput=T expression='!(Y>35503 && Y<35608) && !(Y>30133 && Y<30218) && !(Y>24711 && Y<24827) && !(Y>19314 && Y<19424) && !(Y>13932 && Y<14037) && !(Y<8648) && !(Y>40950) && !(X>39060) && !(X<7650) && !(X>39060)'
    else 
	if ( $instr == 'M1' ) then 
	    echo "MOS1 event cleanup"
	    evselect -w 0 -V 0 table=tempfill.fits filteredset=tempfill2.fits expression='!(Y>20270 && Y<20706 && X>17835 && X<24477) && !(Y>20269 && Y<20377 && X>12162 && X<13117) && !(Y>20271 && Y<33516 && X>17844 && X<17933) && !(Y>33520 && Y<33780 && X>11487 && X<17933) && !(Y>33796 && Y<33744 && X>36830 && X<37518) && !(Y>20429 && Y<20709 && X>24578 && X<28128) && !(Y>20417 && Y<20673 && X>28128 && X<37778) && !(Y>20205 && Y<20350 && X>17017 && X<17817) && !(Y>33675 && Y<33816 && X>31095 && X<37633) && !(Y<20300 && X<11420) && !(Y>33600 && X<11360) && !(Y>33800 && X>37550) && !(Y<20520 && X>37650)' updateexposure=Y writedss=N withfilteredset=Y destruct=Y keepfilteroutput=T
	else 
	    echo "MOS2 event cleanup"
	    evselect -w 0 -V 0 table=tempfill.fits filteredset=tempfill2.fits expression='!(Y>13958 && Y<20400 && X>17693 && X<17906) && !(Y>20385 && Y<20511 && X>17723 && X<30379) && !(Y>20364 && Y<27037 && X>30869 && X<31117) && !(Y>27037 && Y<33730 && X>30904 && X<31130) && !(Y>33730 && Y<40406 && X>30959 && X<31157) && !(Y>33757 && Y<33801 && X>17723 && X<30959) && !(Y>33757 && Y<34542 && X>17745 && X<17795) && !(Y>27049 && Y<27181 && X>30898 && X<40792) && !(Y<13900 && X<17850) && !(Y>40400 && X<17800) && !(Y>40300 && X>31200) && !(Y<13900 && X>31200)' updateexposure=Y writedss=N withfilteredset=Y destruct=Y keepfilteroutput=T
	endif
    endif
	echo 'Move cleaned file'
	if ( -e noattcleaned.fits ) rm -f noattcleaned.fits
	cp tempfill2.fits noattcleaned.fits

	if ( -e noattuncleaned.fits ) rm -f noattuncleaned.fits
	cp tempfill.fits noattuncleaned.fits

	if ( -e uncleaned.fits ) rm -f uncleaned.fits
	cp tempfill.fits uncleaned.fits

	attcalc -V 0 -w 0 eventset=uncleaned.fits attitudelabel=fixed refpointlabel=user withatthkset=N fixedra=$ra_pnt fixeddec=$dec_pnt fixedposangle=$pa_pnt nominalra=$ranom nominaldec=$decnom  
	
    mv -f tempfill2.fits tempfill.fits 

endif
   	   
#-----------------------------------------#
# Use attcalc to recalculate positions
#-----------------------------------------#

echo 'Using position angles '
echo $ra_pnt $dec_pnt $pa_pnt

attcalc -V 0 -w 0 eventset=tempfill.fits attitudelabel=fixed refpointlabel=user withatthkset=N fixedra=$ra_pnt fixeddec=$dec_pnt fixedposangle=$pa_pnt nominalra=$ranom nominaldec=$decnom

attcalc -V 0 -w 0 eventset=$evlist_2 attitudelabel=fixed refpointlabel=user withatthkset=N fixedra=$ra_pnt fixeddec=$dec_pnt fixedposangle=$pa_pnt nominalra=$ranom nominaldec=$decnom

#-----------------------------------------#
# Clear up
#-----------------------------------------#

mv -f tempfill.fits $finalfile
if ($keep == 'Y') mv -f $evlist_2 sourceremoved.fits

rm -f tempfd
#rm -f region.reg
rm -f VERIFY
rm -f $evlist_2
rm -f $myfirstBATCH_file

echo 'Task completed'

exit