combine (version 2.1)

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

###########################################################################
#
# SYNTAX: combine
#
# BRIEFLY: Merge the event files created by frfread.
#
# DESCRIPTION: This routine combines the event files produced by frfread
# DESCRIPTION: to produce unprocessed, unfiltered event files.
# DESCRIPTION: The observation mode of an event file is specified by
# DESCRIPTION: about 70 FITS keywords. Event files which have identical values
# DESCRIPTION: for all of these keywords are merged into a single event
# DESCRIPTION: file for further processing.
# DESCRIPTION: <P>
# DESCRIPTION: Groups of event files with fewer than %gis_events_min (GIS)
# DESCRIPTION: or %sis_events_min (SIS) events are not combined. This is
# DESCRIPTION: because such small files are generally created when ASCA
# DESCRIPTION: is changing modes, and these trasitional mode data should not
# DESCRIPTION: be used for scientific analysis.
# DESCRIPTION: <P>
# DESCRIPTION: After a group is combined, the raw event files in that group
# DESCRIPTION: are deleted. After all groups are combined a number of
# DESCRIPTION: transitional mode "scraps" are left over. Again, these files
# DESCRIPTION: do not contain useful information, so they are not distributed.
# DESCRIPTION: Nevertheless the scraps are tar-ed together and later
# DESCRIPTION: saved as a trend product in case they are needed for
# DESCRIPTION: diagnostic information.
# DESCRIPTION: <P>
# DESCRIPTION: The ORIGMODE FITS keyword is added to each combined event file
# DESCRIPTION: At this stage the value of this keyword is the same as the value
# DESCRIPTION: of the DATAMODE keyword. When SIS FAINT mode files are converted
# DESCRIPTION: to BRIGHT and BRIGHT2 mode, the ORIGMODE keyword can be used
# DESCRIPTION: to distinguish on-board bright mode files from converted
# DESCRIPTION: bright mode files.
# DESCRIPTION: <P>
# DESCRIPTION: Also the TLMIN1 and TLMAX1 keywords are added to each combined
# DESCRIPTION: event file. These give the minimum and maximum values of the
# DESCRIPTION: TIME column and may be needed by downstream software.
#
# VERSION: 2.1
#
# HISTORY: 0.0 -> 1.0 8/23/96
# HISTORY: Now only deletes the raw files which
# HISTORY: are used to make unfiltered files. The remaining "scraps"
# HISTORY: are tarred together to make the scraps file, which is kept as
# HISTORY: a trend product.
# HISTORY:
# HISTORY: 1.0 -> 1.1 9/20/96
# HISTORY: Fixed a bug which used i for the index in two nested loops.
# HISTORY: The result was that only the first telemetry file's raw files
# HISTORY: were being merged. This was fixed by changing the outermost
# HISTORY: loop's index from i to telnum.
# HISTORY:
# HISTORY: 1.1 -> 1.2 9/25/96
# HISTORY: The last fix was still not good enough, since the index i was being
# HISTORY: reset each time the script went to a new telemetry file. The result
# HISTORY: was that unscreened files were getting overwritten.
# HISTORY: Fixed things by saving the running count for each instrument
# HISTORY: individually.
# HISTORY:
# HISTORY: 1.2 -> 1.3 10/8/96
# HISTORY: Modified to use the cmerge STOOL to merge event files instead
# HISTORY: of the merge utility script.
# HISTORY:
# HISTORY: 1.3 -> 1.4 12/10/96
# HISTORY: Corrected a bug where the call to the merge utility was never
# HISTORY: deleted. cmerge was called immediately afterwards, so the cmerge
# HISTORY: output would clobber the merge output, but it must have been
# HISTORY: slowing things down considerably.
# HISTORY:
# HISTORY: 1.4 -> 1.5 1/23/97
# HISTORY: Added a check to see if there are any scraps files to tar before
# HISTORY: trying to tar them. Also fixed bad syntax in the exception call
# HISTORY: for the tar.
# HISTORY: 2/5/97 Added a log entry giving the file names of the merged
# HISTORY: event files.
# HISTORY:
# HISTORY: 1.5 -> 1.6 3/18/97
# HISTORY: Now adds ORIGMODE keyword to all eventfiles. This is to later
# HISTORY: distinguish on-board- and ground-converted BRIGHT mode files.
# HISTORY:
# HISTORY: 1.6 -> 1.7 6/25/97
# HISTORY: Now adds the TLMIN1 and TLMAX1 values to the first extension of the
# HISTORY: combined event files, as requested by Keith Arnaud.
# HISTORY:
# HISTORY: 1.7 -> 1.8 6/26/97
# HISTORY: Commented out the last change, since it screwed up the faint
# HISTORY: FTOOL.
# HISTORY:
# HISTORY: 1.8 -> 1.9 9/25/97
# HISTORY: Added some sanity checks for bad event files.
# HISTORY: Made more quiet when there are no raw files for a given instrument.
# HISTORY:
# HISTORY: 1.9 -> 1.10 10/9/97
# HISTORY: Added a check so that RAWXBINS=RAWYBINS is only checked for PHA mode.
# HISTORY: Also resolved a variable name collision for "telem" that was created
# HISTORY: by adding the check for "falling through the cracks" of the attitude
# HISTORY: file.
# HISTORY:
# HISTORY: 1.10 -> 2.0 10/16/97
# HISTORY: Added code to check the OBJECT keywords in the unfiltered and
# HISTORY: housekeeping files. There is apparently a bug in frfread3.032
# HISTORY: which doesn't handle commas in object names correctly.
# HISTORY: If the keywords differ from the value in the jobpar, then
# HISTORY: all the keywords are fixed.
# HISTORY: 12/12/97
# HISTORY: Fixed a bug which was expecting DATAMODE="PHA" instead of "PH"
# HISTORY: This caused GIS files to not be checked for valid RAWXBINS
# HISTORY: and RAWYBINS keywords.
# HISTORY: Fixed a bug that wasn't catching the "fall through the cracks
# HISTORY: of the attitude file" problem.
# HISTORY: Added a check for if SIS SPTR keywords exist.
# HISTORY:
# HISTORY: 2.0 -> 2.1 4/24/98
# HISTORY: Added code to handle both compressed and gzipped telemetry files.
# HISTORY: Now do not combine GIS PCAL mode files. These are rare and can not
# HISTORY: be analyzed using the current software. It is also likely that
# HISTORY: they have not been well tested in frfread and may contain bugs.
# HISTORY: Their raw files will be saved with the
# HISTORY: scraps.
#
############################################################################
#DEBUG=1

$UTIL/milestone "Merging event files from frfread"

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

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


################################################################
# Loop through the parent telemetry files so that
# event files from different telemetry files will not be merged
################################################################

##################################
# Initialize running file indices
##################################
sis0index=1
sis1index=1
gis2index=1
gis3index=1
typeset -i i

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

    ####################################
    # Get base filename for event files
    ####################################
    number="$(echo $telnum | nawk '{printf ("%03d", $0)}')"
    telem=$($UTIL/read_parfile $JOBPAR frffile_$number)
    telem=${telem##*/}
    telem=${telem%.Z}
    telem=${telem%.gz}
    telem_base=$(echo $telem | tr '.' '_')

    ###################
    # debugging output
    ###################
    if [ -n "$DEBUG" ]; then
        echo ${0##*/}: telnum=$telnum
        echo ${0##*/}: telem=$telem
    fi


    for inst in g2 g3 s0 s1; do

        case "$inst" in
        g2)
            caps=G2
            display=GIS2
            min_events=$($UTIL/read_parfile $PARFILE gis_events_min );;
        g3)
            caps=G3
            display=GIS3
            min_events=$($UTIL/read_parfile $PARFILE gis_events_min );;
        s0)
            caps=S0
            display=SIS0
            min_events=$($UTIL/read_parfile $PARFILE sis_events_min );;
        s1)
            caps=S1
            display=SIS1
            min_events=$($UTIL/read_parfile $PARFILE sis_events_min );;
        esac

        #########################
        # Sort files by keywords
        #########################
        $UTIL/log_entry "Collecting $display event files by mode"

        case $inst in
        g?)
            ############
            # GIS files
            ############
            ls ${telem_base}${caps}*7[01][HML].fits 2>/dev/null | \
                             $STOOL/gissortcode >code.tmp 2>stderr.log
            $UTIL/sort_test gissortcode $0 1 stderr.log

            sort code.tmp | $STOOL/gissortsplit >split.tmp 2>stderr.log
            $UTIL/sort_test gissortsplit $0 1 stderr.log;;
        s0)
            #############
            # SIS0 files
            #############
            ls ${telem_base}${caps}*[HML].fits 2>/dev/null | \
                             $STOOL/sis0sortcode >code.tmp 2>stderr.log
            $UTIL/sort_test sis0sortcode $0 1 stderr.log

            sort code.tmp | $STOOL/sis0sortsplit >split.tmp 2>stderr.log
            $UTIL/sort_test sis0sortsplit $0 1 stderr.log;;
        s1)
            #############
            # SIS1 files
            #############
            ls ${telem_base}${caps}*[HML].fits 2>/dev/null | \
                             $STOOL/sis1sortcode >code.tmp 2>stderr.log
            $UTIL/sort_test sis1sortcode $0 1 stderr.log

            sort code.tmp | $STOOL/sis1sortsplit >split.tmp 2>stderr.log
            $UTIL/sort_test sis1sortsplit $0 1 stderr.log;;
        esac

        rm -f code.tmp

        sort -r split.tmp >listlist.tmp
        rm -f split.tmp

        #########################################################
        # Set the index to the correct value for that instrument
        #########################################################
        case "$inst" in
        s0) i=$sis0index;;
        s1) i=$sis1index;;
        g2) i=$gis2index;;
        g3) i=$gis3index;;
        esac

        maxevents=$(awk '{print $1}' listlist.tmp | sort -n | tail -1 )

        for list in  $(awk '{print $2}' listlist.tmp); do

            nevents="$(grep $list listlist.tmp | awk '{print $1}')"

            #######################################################
            # Check if this is a valid group of event files, or if
            # it's just scraps
            #######################################################
            keepit="yes"

            ###################
            # Debugging output
            ###################
            if [ -n "$DEBUG" ]; then
                echo "${0##*/}: maxevents=$maxevents"
                echo "${0##*/}: nevents=$nevents"
                echo "${0##*/}: wc -l <$list=$(wc -l <$list)"
                echo "contents of $list"
                cat $list
                echo "End of list"
            fi

            if [ $nevents -lt $min_events ]; then
                keepit="no"
            else
                if [ $nevents -lt $maxevents -a \
                     $(wc -l <$list) -eq 1 ]; then
                    ######################################################
                    # The proposed combined file has fewer than the
                    # maximum number of events, and it would be composed
                    # of only one raw file.
                    # Make some extra sanity checks to make sure the file
                    # isn't weird.
                    ######################################################
                    file=$(cat $list )
                    $UTIL/log_entry "Checking $file"

                    datamode=$($UTIL/read_fits_keyword $file DATAMODE)

                    #####################
                    # Check GIS keywords
                    #####################
                    case "$datamode" in
                    PH)
                        ########################
                        # RAWXBINS = RAWYBINS ?
                        ########################
                        rawxbins=$($UTIL/read_fits_keyword $file RAWXBINS)
                        rawybins=$($UTIL/read_fits_keyword $file RAWYBINS)

                        if [ $rawxbins -ne $rawybins ]; then
                            $UTIL/exception $0 1 \
                            "RAWXBINS=$rawxbins and RAWYBINS=$rawybins in $file"
                            keepit="no"
                        fi
                        ;;

                    BRIGHT|BRIGHT2)
                        #######################################################
                        # check SIS keywords
                        # specifically, check that the SPTR keywords
                        # exist for all active chips.
                        # The split threshold for a chip is output to the
                        # telemetry once every 4 superframes, so short files
                        # at the beginning of the observation may legitimately
                        # be missing these keywords. We could just add SPLT
                        # keyword, but many other keywords will also be
                        # missing which are not easy to replace since their
                        # values may change within the observation. Also
                        # such files are necessarily small, so it's no great
                        # loss
                        #######################################################
                        ccdlist=$($UTIL/read_fits_keyword $file[0] \
                                                          S${inst#s}CCDLST )
                        for chip in $ccdlist; do

                            key=S${inst#s}_SPTR${chip}

                            exist=$($UTIL/keyword_exists $file[0] $key )
                            if [ "$exist" != "yes" ]; then
                                ############################################
                                # keyword is missing
                                ############################################
                                $UTIL/exception $0 1 \
                                     "Deleting $file with missing $key keyword"
                                keepit=no
                                break
                            fi

                        done
                        ;;
                    esac

                    ################################################
                    # Fall through the cracks of the attitude file?
                    ################################################
                    tstart=$($UTIL/read_fits_keyword $file TSTART )
                    tstop=$( $UTIL/read_fits_keyword $file TSTOP  )

                    duration=$($STOOL/equals $tstop - $tstart )
                    if [ $($STOOL/floatcmp $duration 100. ) = "lt" ]; then
                        ###############################################
                        # Event file duration shorter than 100 kilosec
                        # Better check this one
                        ###############################################
                        att=fa${telem#ft}
                        rm -f time.tmp
                        $UTIL/dump_table $att TIME >time.tmp

                        rows=$(nawk '$1 >= time ' time=$tstart time.tmp | \
                               nawk '$1 <= time ' time=$tstop | \
                               wc -l)

                        rm -f time.tmp

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

                        if [ $rows -eq 0 ]; then
                             $UTIL/exception $0 1 \
                                            "No attitude rows for $file"
                             keepit="no"
                        fi
                    fi
                fi
            fi

            #####################################
            # Combine the files if they are good
            #####################################
            if [ "$keepit" = "yes" ]; then
                ######################################
                # This is a good group, so combine it
                ######################################
                index="$(echo $i | nawk '{printf ("%03d", $0)}')"

                mode=${list#?????}
                mode=${mode%?.prelist}

                bitrate=${list#???????}
                bitrate=${bitrate%.prelist}

                unfiltered=$($UTIL/generate_filename unfiltered \
                             $inst $index $mode $bitrate )

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

                ##############
                # Merge files
                ##############
                rm -f $unfiltered
                $UTIL/log_entry "Creating $unfiltered"
                $STOOL/cmerge $unfiltered @$list >stdout.log 2>stderr.log
                $UTIL/stool_test cmerge $? $0 2 stdout.log stderr.log

                #######################################################
                # Add ORIGMODE keyword with same value as DATAMODE
                # This is to later distinguish ground converted BRIGHT
                # mode from on-board BRIGHT mode
                #######################################################
                datamode=$($UTIL/read_fits_keyword $unfiltered DATAMODE )

                $UTIL/add_fits_keyword $unfiltered[0] ORIGMODE $datamode \
                                       "DATAMODE before any conversion"

                $UTIL/add_fits_keyword $unfiltered[1] ORIGMODE $datamode \
                                       "DATAMODE before any conversion"

                $UTIL/add_fits_keyword $unfiltered[2] ORIGMODE $datamode \
                                       "DATAMODE before any conversion"

                ###################################################
                # Add TLMIN and TLMAX keywords for the time column
                ###################################################
                tstart=$($UTIL/read_fits_keyword $unfiltered TSTART )
                tstop=$( $UTIL/read_fits_keyword $unfiltered TSTOP  )

                #######################
                # Remove the raw files
                #######################
                for file in $(cat $list ); do

                     rm -f $file

                done

                let i=$i+1

            else
                ###################################
                # Not enough events to bother with
                ###################################
                $UTIL/log_entry \
                    "Ignoring the following files containing $nevents events"
                $UTIL/file_to_log $list

            fi

            rm -f $list

        done #loop over unscreened files

        rm -f listlist.tmp

        ##################################
        # Remember the running file index
        ##################################
        case "$inst" in
        s0) sis0index=$i;;
        s1) sis1index=$i;;
        g2) gis2index=$i;;
        g3) gis3index=$i;;
        esac

    done #loop over instruments

    let telnum=$telnum+1

done # loop over telemetry files

##########################
# Tar together the scraps
##########################
$UTIL/log_entry "Tar-ing together the leftover raw files"
scraps=$($UTIL/generate_filename scraps )

list=$(ls ft*[HML].fits 2>/dev/null)
if [ -n "$list" ]; then
    tar -cvf $scraps ft*[HML].fits >stdout.log 2>stderr.log
    code=$?

    $UTIL/file_to_log stdout.log
    if [ $code -ne 0 -o -s stderr.log ]; then

        $UTIL/exception $0 1 "Error in tar, exit code=$code"

        if [ -s stderr.log ]; then
            $UTIL/log_entry "Stderr output from tar:"
            $UTIL/file_to_log stderr.log
        fi
    fi

    rm -f stdout.log
    rm -f stderr.log

    rm -f ft*[HML].fits
fi


####################################################################
# Check that the OBJECT keyword is correct in the event and HK file
# headers. There is a bug in frfread3.032 which has trouble with
# object names containing commas
####################################################################
$UTIL/log_entry "Checking OBJECT keywords in HK and event files"

file=$(ls ft*HK.fits 2>/dev/null | head -1 )
keyword=$($UTIL/read_fits_keyword $file[0] OBJECT )

file=$($UTIL/read_parfile $JOBPAR frffile_001 )
file=${file##*/}
file=${file%.Z}
file=${file%.gz}
object=$($UTIL/read_fits_keyword $file[1] OBJECT)

if [ "$object" != "$keyword" ]; then
    ###############################
    # Need to fix all the keywords
    ###############################
    $UTIL/exception $0 1 "OBJECT keywords $keyword do not match $object"
    $UTIL/log_entry "Fixing OBJECT keywords in HK and unfiltered files"

    #####################
    # Housekeeping files
    #####################
    for file in $(ls ft*HK.fits 2>/dev/null ); do

        for ext in 0 1; do

            $UTIL/add_fits_keyword $file[$ext] OBJECT "$object"

        done

    done

    ##############
    # Event files
    ##############
    list=$($UTIL/generate_filename unfiltered any any any any any )
    list=$(ls $list 2>/dev/null )
    for file in $list; do

        for ext in 0 1 2; do

            $UTIL/add_fits_keyword $file[$ext] OBJECT "$object"

        done

    done

fi


exit 0