linearize (version 1.7)

You can also look at:
#! /usr/bin/ksh

###########################################################################
#
# SYNTAX: linearize
#
# BRIEFLY: Run ascalin to calibrate the unfiltered event files.
#
# DESCRIPTION: This routine runs the <TT>ascalin</TT>
# DESCRIPTION: FTOOL on all of the unfiltered
# DESCRIPTION: event files. This does several things.
# DESCRIPTION: PI values are calculated from PHA (energy) values.
# DESCRIPTION: For the SIS this means normalizing the four chips to a common
# DESCRIPTION: scale. For the GIS, corrections are made for the temperature
# DESCRIPTION: dependence of the gain.
# DESCRIPTION: To do this ascalin uses a gain history file which is
# DESCRIPTION: created by the ISAS program <TT>temp2gain</TT>.
# DESCRIPTION: <P>
# DESCRIPTION: For the GIS, raw X and Y coordinates of events are converted
# DESCRIPTION: to a linear detector X and Y scale.
# DESCRIPTION: For both instruments, the detector coordinates are converted
# DESCRIPTION: to R.A. and Dec. coordinates registered about a common aspect
# DESCRIPTION: point.
#
# VERSION: 1.7
#
# HISTORY: 0.0 -> 1.0 12/19/96
# HISTORY: Now reads RA and dec from the jobpar.
# HISTORY: 12/23/96 Also skip MPC mode event files.
# HISTORY: 12/27/96 Take name of temp2gain executable from the parfile.
# HISTORY: Assume that the temp2gain parfile has the same name and location
# HISTORY: as the executable, but with ".par" appended.
# HISTORY:
# HISTORY: 1.0 -> 1.1 6/13/97
# HISTORY: Now removes unfiltered file if there was an error from ascalin.
# HISTORY: If a file can't be calibrated, it can't be used, and it can
# HISTORY: seriously screw up downstream processing, so this is the best way
# HISTORY: to recover, although it does make debugging harder.
# HISTORY:
# HISTORY: 1.1 -> 1.2 7/10/97
# HISTORY: Now just renames files for which ascascreen has an error.
# HISTORY: This has the same effect as removing the file, but makes
# HISTORY: debugging easier.
# HISTORY:
# HISTORY: 1.2 -> 1.3 9/2/97
# HISTORY: Fixed bug when moving ascalin-failed files to "ad*.tmp" files
# HISTORY: in order to remove them. The bug resulted in creating ".tmp"
# HISTORY: files which remained hidden until ascajb choked on them, since
# HISTORY: they weren't cataloged.
# HISTORY: 9/16/97 Now also does a better job handling "out of time order"
# HISTORY: errors from ascalin. Before, this error wasn't caught, and
# HISTORY: the event file was deleted. Now this sort of error is caught
# HISTORY: and not reported, and the event file is sorted in time.
# HISTORY:
# HISTORY: 1.3 -> 1.4 9/25/97
# HISTORY: Added a catch to make things go smoother if there are no unfiltered
# HISTORY: event files. Also, add .unf files which fail ascalin to the
# HISTORY: end of the scraps file.
# HISTORY:
# HISTORY: 1.4 -> 1.5 12/16/97
# HISTORY: Added code to check each unfiltered file before it is run through
# HISTORY: ascalin to see if the file is covered by the attitude file.
# HISTORY: If it isn't, the uncovered parts are clipped off and saved
# HISTORY: among the scraps. This replaces a previous method of dealing
# HISTORY: with this, which trimmed the telemetry file to have the same
# HISTORY: start and stop times as the attitude file. This did not always
# HISTORY: work, though, since event times can extend beyond the last time
# HISTORY: in the telemetry file. That was also bad because it modified
# HISTORY: the original telemetry and it threw away some potentially
# HISTORY: usable photons.
# HISTORY:
# HISTORY: 1.5 -> 1.6 4/08/98
# HISTORY: Changed "column" to "columns" in fmemsort to accomodate FTOOLS 4.1
# HISTORY: Also, now code handles both gzipped and compressed telemetry and
# HISTORY: attitude files.
# HISTORY: Added code to calculate Eric Gothelf's attitude corrections
# HISTORY: and store the results in the job par file.
# HISTORY:
# HISTORY: 1.6 -> 1.7 1998-08-07
# HISTORY: Modified to use temp2gain 4.3
# HISTORY: Also deleted reference to ISASBIN
#
############################################################################
#DEBUG=1

$UTIL/milestone "Creating GIS gain history file"

##############################################
# Read number of telemetry files from job par
##############################################
nfiles=$($UTIL/read_parfile $JOBPAR frfcnt)
if [ -z "$nfiles" ]; then
    nfiles=0
fi

if [ $nfiles -eq 0 ]; then
    $UTIL/exception $0 3 "Can't get non-zero number of frf files from $JOBPAR"
fi

###############################
# Loop through telemetry files
###############################
caldir=$($UTIL/read_parfile $PARFILE caldir)
calfile=$($UTIL/read_parfile $PARFILE gainhist)
leaptable=$($UTIL/read_parfile $PARFILE leaptable)

$UTIL/fetch_file $caldir $calfile
$UTIL/fetch_file $caldir $leaptable

SYSPFILES_save="$SYSPFILES"
PFILES_save="$PFILES"


typeset -i i=1
while [ $i -le $nfiles ]; do

    number="$(echo $i | nawk '{printf ("%03d", $0)}')"
    full=$($UTIL/read_parfile $JOBPAR frffile_$number)
    file=${full##*/}
    file=${file%.Z}
    file=${file%.gz}

    ghf=$(echo $file | tr '.' '_').ghf

    $UTIL/log_entry "Creating gain history file $ghf"

    temp2gain=$($UTIL/read_parfile $PARFILE temp2gain)

    $UTIL/setup_parfile $temp2gain.par frffile="$file"        \
                                       outfile="$ghf"         \
                                       tblfile="${file}.tbl"  \
                                       leapsec="./$leaptable" \
                                       type="0"               \
                                       steptime="600"         \
                                       histfile=$calfile

    rm -f $ghf

    ########################################################################
    # temporarily reset the parfile environment variables to the temp2gain
    # directory, and then set them back before any FTOOLS are read
    #######################################################################
    export SYSPFILES="${temp2gain%/*}"
    export PFILES="$LOCPFILES;$SYSPFILES"

    $temp2gain >stdout.log 2>stderr.log
    code=$?

    export SYSPFILES="SYSPFILES_save"
    export PFILES="$PFILES_save"


    errorlines=$(grep "Error|Warning" stdout.log | wc -l)

    if [ $code -ne 0 -o $errorlines -gt 0 ]; then
        $UTIL/exception $0 2 \
                        "Error from ISAS program temp2gain. Exit code=$code"
    fi

    $UTIL/log_entry "Output from ISAS program temp2gain"
    $UTIL/file_to_log stderr.log
    if [ -s stdout.log ]; then
        $UTIL/log_entry "Stdout output from ISAS program temp2gain"
        $UTIL/file_to_log stdout.log
    fi

    rm -f stdout.log stderr.log
    rm -f temp2gain.par
    rm -f ${file}.tbl

    ####################################
    # Check if ghf file is time ordered
    ####################################
    $UTIL/check_order $ghf[1] CAL_START
    $UTIL/check_order $ghf[1] CAL_STOP

    #################################
    # Tell the ghf files their names
    #################################
    $UTIL/add_fits_keyword $ghf[0] FNAME $ghf
    $UTIL/add_fits_keyword $ghf[1] FNAME $ghf

    let i=$i+1

done

########################################
########################################
#
# calculate attitude adjustments
#
########################################
########################################
i=1
while [ $i -le $nfiles ]; do

    number="$(echo $i | nawk '{printf ("%03d", $0)}')"
    full=$($UTIL/read_parfile $JOBPAR frffile_$number)
    telem=${full##*/}
    telem=${telem%.Z}
    telem=${telem%.gz}

    chk=$(echo $telem | tr '.' '_')CMHK.fits
    att=$(echo $telem | tr 't' 'a')

    tstart=$($UTIL/read_fits_keyword $att[1] MTIME0)
    tstop=$( $UTIL/read_fits_keyword $att[1] MTIME1)


    $UTIL/log_entry "Calculating attitude correction from $chk"

    if [ ! -a "$file" ]; then
         $UTIL/exception $0 2 "Common HK file $chk does not exist"
         let i=$i+1
         continue
    fi

    $UTIL/attitude_adjustment $chk $tstart $tstop | \
    read mean sigma delx dely delra deldec

    $UTIL/add_parameter $JOBPAR bp4mean_$number "$mean" "r" \
                        "Mean of BP4 HK parameter"

    $UTIL/add_parameter $JOBPAR bp4sigma_$number "$sigma" "r" \
                        "Std Dev of BP4 HK parameter"


    $UTIL/add_parameter $JOBPAR delx_$number "$delx" "r" \
                        "Pointing Error in X (mm)"
    $UTIL/add_parameter $JOBPAR dely_$number "$dely" "r" \
                        "Pointing Error in Y (mm)"

    $UTIL/add_parameter $JOBPAR delra_$number "$delra" "r" \
                        "Pointing Error in R. A. (arcmin)"

    $UTIL/add_parameter $JOBPAR deldec_$number "$deldec" "r" \
                        "Pointing Error in Dec. (arcmin)"

    let i=$i+1

done




###################################
###################################
##
##  ASCALIN
##
###################################
###################################
$UTIL/milestone "Running ASCALIN on unfiltered event files"

frftimetoll=$($UTIL/read_parfile $PARFILE frftimetoll)

##############################################
# Read number of telemetry files from job par
##############################################
nfiles=$($UTIL/read_parfile $JOBPAR frfcnt)
if [ -z "$nfiles" ]; then
    nfiles=0
fi

if [ $nfiles -eq 0 ]; then
    $UTIL/exception $0 3 "Can't get non-zero number of frf files from $JOBPAR"
fi

#############################
# Calibration file directory
#############################
caldir=$($UTIL/read_parfile $PARFILE caldir)

##########################################
# Loop through the unfiltered event files
##########################################
list=$($UTIL/generate_filename unfiltered any any any any)
list=$(ls $list 2>/dev/null)
for unfiltered in ${list} ; do

    inst=$(   $UTIL/parse_filename inst    $unfiltered)
    index=$(  $UTIL/parse_filename index   $unfiltered)
    mode=$(   $UTIL/parse_filename mode    $unfiltered)
    bitrate=$($UTIL/parse_filename bitrate $unfiltered)

    ###############################
    # Skip FAST and MPC mode files
    ###############################
    if [ "$mode" = "03" -o "$mode" = "71" ]; then
        continue
    fi

    ###########################################################
    # Make sure the event file is covered by the attitude file
    ###########################################################
    $UTIL/log_entry "Checking if $unfiltered is covered by attitude file"

    tstart=$($UTIL/read_fits_keyword $unfiltered[0] TSTART)
    tstop=$( $UTIL/read_fits_keyword $unfiltered[0] TSTOP )
    telem=$( $UTIL/read_fits_keyword $unfiltered[0] TLM_FILE)

    att=fa${telem#ft}
    att_start=$($UTIL/read_fits_keyword $att[1] MTIME0)
    att_stop=$( $UTIL/read_fits_keyword $att[1] MTIME1)

    att_start=$($STOOL/equals "$att_start - $frftimetoll")
    att_stop=$( $STOOL/equals "$att_stop  + $frftimetoll")

    ###################
    # Debugging output
    ###################
    if [ -n "$DEBUG" ]; then
        echo ${0##*/}: tstart=$tstart
        echo ${0##*/}: tstop=$tstop
        echo ${0##*/}: att=$att
        echo ${0##*/}: att_start=$att_start
        echo ${0##*/}: att_stop=$att_stop
    fi

    if [ $($STOOL/floatcmp $tstart $att_start) = "lt" -o \
         $($STOOL/floatcmp $tstop  $att_stop ) = "gt" ]; then
        #########################################
        # Just skip over FAINT mode files
        # They will remain, but be un-ascalinned
        #########################################
        if [ "$mode" = "01" ]; then
            $UTIL/exception $0 1 \
                        "$unfiltered not ascalined since not covered by $att"
            continue
        fi

        #################################
        # Create some ascii time filters
        #################################
        rm -f covered.tmp
        echo "$att_start $att_stop" >covered.tmp

        rm -f uncovered.tmp
        echo "$tstart $att_start" >uncovered.tmp
        echo "$att_stop $tstop"  >>uncovered.tmp

        uncovered=${unfiltered}_no_att

        ###############################################
        # Use the extractor to clip away the bad parts
        ###############################################
        $UTIL/log_entry "Clipping off the uncovered parts of $unfiltered"
        rm -f events.tmp
        $UTIL/setup_parfile $FTOOL/extractor.par \
                            filename=$unfiltered \
                            eventsout=events.tmp  \
                            regionfile="NONE" \
                            qdpfile=NONE \
                            fitsbinlc="NONE" \
                            unbinlc=NONE \
                            phafile=NONE \
                            imgfile=NONE \
                            timefile="covered.tmp" \
                            xcolf=DETX \
                            ycolf=DETY \
                            xcolh=DETX \
                            ycolh=DETY \
                            xfkey=RAWXBINS \
                            yfkey=RAWYBINS \
                            xhkey=RAWXBINS \
                            yhkey=RAWYBINS \
                            phamax=NAXIS \
                            specbin=1 \
                            binh=1 \
                            binf=1 \
                            binlc=1.6000000E+01 \
                            tcol=TIME \
                            ecol=PI \
                            gti="STDGTI" \
                            events="EVENTS" \
                            gtitxt=NONE \
                            gtinam="STDGTI"

        $FTOOL/extractor >stdout.log 2>stderr.log
        $UTIL/ftool_test extractor $? $0 2 stdout.log stderr.log

        ####################################################
        # Now extract the bad parts to save with the scraps
        ####################################################
        $UTIL/log_entry "Extracting the uncovered parts to $uncovered"

        rm -f $uncovered
        $UTIL/setup_parfile $FTOOL/extractor.par \
                            filename=$unfiltered \
                            eventsout=$uncovered  \
                            regionfile="NONE" \
                            qdpfile=NONE \
                            fitsbinlc="NONE" \
                            unbinlc=NONE \
                            phafile=NONE \
                            imgfile=NONE \
                            timefile="uncovered.tmp" \
                            xcolf=DETX \
                            ycolf=DETY \
                            xcolh=DETX \
                            ycolh=DETY \
                            xfkey=NAXIS \
                            yfkey=NAXIS \
                            xhkey=NAXIS \
                            yhkey=NAXIS \
                            phamax=NAXIS \
                            specbin=1 \
                            binh=1 \
                            binf=1 \
                            binlc=1.6000000E+01 \
                            tcol=TIME \
                            ecol=PI \
                            gti="STDGTI" \
                            events="EVENTS" \
                            gtitxt=NONE \
                            gtinam="STDGTI"

        $FTOOL/extractor >stdout.log 2>stderr.log
        $UTIL/ftool_test extractor $? $0 2 stdout.log stderr.log

        ###########################################
        # Save the uncovered part among the scraps
        ###########################################
        $UTIL/log_entry "Adding $uncovered to the scraps"
        scraps=$($UTIL/generate_filename scraps)
        tar -uvf $scraps $uncovered >stdout.log
        $UTIL/file_to_log stdout.log
        rm -f stdout.log

        ##################################################
        # Replace the unfiltered file by the covered part
        ##################################################
        nevents=$($UTIL/read_fits_keyword events.tmp[1] NAXIS2)
        if [ $nevents -gt 0 ]; then
            mv -f events.tmp $unfiltered
        else
            $UTIL/log_entry "No events in $unfiltered were covered by $att"
            continue
        fi

        ################
        # Error message
        ################
        ontime=$($UTIL/read_fits_keyword $uncovered[0] ONTIME)
        ontime=$($STOOL/equals "int( $ontime +.5)")
        $UTIL/exception $0 1 \
                "Removed $ontime seconds from $unfiltered not covered by $att"

        rm -f covered.tmp
        rm -f uncovered.tmp
        rm -f $uncovered
    fi

    ############
    # Log entry
    ############
    $UTIL/log_entry "Running ascalin on $unfiltered"

    teldef=$($UTIL/generate_filename teldef $inst $index $mode $bitrate)
    telemetry=$($UTIL/generate_filename telemetry $inst $index $mode $bitrate)

    case "$inst" in
    g*)
        ######
        # GIS
        ######
        tempofile="$(echo ${telemetry} | tr '.' '_').ghf";;
    s*)
        ######
        # SIS
        ######
        tempofile=$($UTIL/read_parfile $PARFILE sispical)
        $UTIL/fetch_file $caldir $tempofile;;
    esac

    #####################################################
    # Find attitude file corresponding to telemetry file
    #####################################################
    file=""
    typeset -i i=1;
    while [ $i -le $nfiles -a "$file" != "$telemetry" ]; do

        number="$(echo $i | nawk '{printf ("%03d", $0)}')"

        file=$($UTIL/read_parfile $JOBPAR frffile_$number)
        file=${file##*/}
        file=${file%.Z}
        file=${file%.gz}


        ###################
        # Debugging output
        ###################
        if [ -n "$DEBUG" ]; then
            echo ${0##*/}: i=$i
            echo ${0##*/}: file=$file
        fi

        let i=$i+1

    done

    if [ "$file" = "$telemetry" ]; then
        attitude=$($UTIL/read_parfile $JOBPAR attfile_$number)
        attitude=${attitude##*/}
        attitude=${attitude%.Z}
        attitude=${attitude%.gz}

    else
        $UTIL/exception $0 3 "No attitude file corresponding to $telemetry"
        attitude=none
    fi

    #############
    # RA and Dec
    #############
    ra=$( $UTIL/read_parfile $JOBPAR ra)
    dec=$($UTIL/read_parfile $JOBPAR dec)

    ###################
    # Debugging output
    ###################
    if [ -n "$DEBUG" ]; then
        echo ${0##*/}: unfiltered=$unfiltered
        echo ${0##*/}: index=$index
        echo ${0##*/}: mode=$mode
        echo ${0##*/}: bitrate=$bitrate
        echo ${0##*/}: teldef=$teldef
        echo ${0##*/}: telemetry=$telemetry
        echo ${0##*/}: attitude=$attitude
    fi

    $UTIL/fetch_file $caldir $teldef

    $UTIL/setup_parfile $FTOOL/ascalin.par calfile=$teldef      \
                                           tempofile=$tempofile \
                                           attitude=$attitude   \
                                           pointing=USER        \
                                           datafile=$unfiltered \
                                           ranom=$ra            \
                                           decnom=$dec          \
                                           history=no           \
                                           verbose=no

    $FTOOL/ascalin >stdout.log 2>stderr.log
    code=$?
    grep -v "Detected gap > 15min in attitude file:" stderr.log | \
    grep -v "out of time order" >filt.log

    dum=$(grep "out of time order" stderr.log | wc -w)
    if [ $dum -ne 0 ]; then
        ################################################
        # Sort the event file since it was out of order
        ################################################
        $UTIL/log_entry "Sorting $unfiltered"
        rm -f out.tmp
        $UTIL/setup_parfile  $FTOOL/fmemsort.par infile=$unfiltered  \
                                                 outfile=out.tmp     \
                                                 columns=TIME        \
                                                 method=insert       \
                                                 ascend=yes          \
                                                 copyall=yes         \
                                                 history=yes
        $FTOOL/fmemsort >stdout_sort.log 2>stdout_sort.log
        $UTIL/ftool_test fmemsort $? $0 2 stdout_sort.log stderr_sort.log
    fi

    $UTIL/ftool_test ascalin $code $0 2 stdout.log stderr.log filt.log

    if [ $? -ne 0 ]; then
        ###############################################
        # Remove unfiltered file if there was an error
        ###############################################
        $UTIL/log_entry "Removing $unfiltered since it can't be calibrated"
        scraps=$($UTIL/generate_filename scraps)
        tar -uvf $scraps $unfiltered >stdout.log
        $UTIL/file_to_log stdout.log
        rm -f stdout.log
        rm -f  $unfiltered
    fi

done


exit 0