Subs::Filter (version 2.0)


package Subs::Filter;
##############################################################################
#
# DESCRIPTION: Filter the event data. First calculate attitude and
# DESCRIPTION: orbit related filtering quantities using the prefilter
# DESCRIPTION: tool. Then create the filter file by collating various
# DESCRIPTION: sources of HK with the prefilter output.
# DESCRIPTION: Next apply filtering criteria - spacial, temporal, or
# DESCRIPTION: event selection - to the event data.
# DESCRIPTION:
#
# HISTORY:
# HISTORY: $Log: Filter.pm,v $
# HISTORY: Revision 1.29  2005/04/06 15:45:03  apsop
# HISTORY: Change to using CALDB for cal parameters.
# HISTORY:
# HISTORY: Revision 1.28  2005/03/25 19:55:03  apsop
# HISTORY: Fix bug in merging of HK data into mkf file.
# HISTORY:
# HISTORY: Revision 1.27  2005/03/15 19:57:56  apsop
# HISTORY: Fix bug in call for getting eng hk file names.  Write DATE-* keywords to mkf primary header.
# HISTORY:
# HISTORY: Revision 1.27  2005/03/15 19:03:26  apsop
# HISTORY: Check for presence of xrt no position message when making tdrss cat file.
# HISTORY:
# HISTORY: Revision 1.26  2005/02/08 18:18:38  apsop
# HISTORY: Move calculation of attorb file to Filter class. Split making of attitude flag columns in to two parts.
# HISTORY:
# HISTORY: Revision 1.25  2005/01/24 18:58:11  apsop
# HISTORY:
# HISTORY: Add support for SAFEHOLD column to makefilter file.
# HISTORY:
# HISTORY: Revision 1.24  2005/01/12 17:21:02  apsop
# HISTORY: Incorporate removal of _PNT suffix from output columns of prefilter
# HISTORY:
# HISTORY: Revision 1.23  2004/11/19 21:46:48  apsop
# HISTORY: New version of xrt2fits.
# HISTORY:
# HISTORY: Revision 1.22  2004/10/12 16:19:49  apsop
# HISTORY: Turn on data filtering to produce cleaned event lists
# HISTORY:
# HISTORY: Revision 1.21  2004/08/30 13:21:51  apsop
# HISTORY: Bobs new version of module, but comment out actual screening for now.
# HISTORY:
# HISTORY: Revision 1.20  2004/08/27 18:38:13  apsop
# HISTORY: Modified to use Ed's code to determine the filtering expressions, but
# HISTORY: the instrument tools for the actual filtering.
# HISTORY:
# HISTORY: Revision 1.19  2004/05/06 20:02:34  dah
# HISTORY: Add version number back into the header comments.
# HISTORY:
# HISTORY: Revision 1.18  2004/04/16 20:21:18  dah
# HISTORY: Begin using embedded history records
# HISTORY:
#
# VERSION: 2.0
#
#
##############################################################################


use Subs::Sub;

@ISA = ("Subs::Sub");
use strict;

sub new {
    my $proto=shift;
    my $self=$proto->SUPER::new();

    $self->{DESCRIPTION}="Cleaning and filtering the event files";

    return $self;
}

##################
# METHODS:
##################

sub body {
    my $self=shift;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();

    #########################
    # create the filter file
    #########################
    $self->make_filter_file();

    ##################################################
    # plot interesting quantities in the filter file
    ##################################################
  #  $self->make_plots();

    ##################################################
    # screen out bad event data
    ##################################################
    $self->filter_data();


} # end of body method

#############################################################################
# assemble everything together into a single filter file
#############################################################################
sub make_filter_file {
    my $self=shift;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();

    #############################################################
    # get the name of the filter file
    #############################################################
    my $filter = $filename->get("filter", "s");
    unlink $filter;
    
    $log->entry("Creating filter file $filter");

    ####################################
    # create the configuration file
    ####################################
    my $config = $self->temp_file("makefilter_config");
    unlink $config;
    open CONFIG, ">$config";

    #########################################################
    # collect attitude/orbit derived information
    #########################################################
    my $attorb = $filename->get("attorb", "s");
    if( -f $attorb) {
        print CONFIG "ELV $attorb PREFILTER D D % ".
                     "angle between Z axis and earth limb\n";

        print CONFIG "BR_EARTH $attorb PREFILTER D D % ".
                     "Angular distance from the sunlit earth\n";
                     
        print CONFIG "FOV_FLAG $attorb PREFILTER D D % ".
                     "0=sky; 1=dark earth; 2=bright earth\n";

        print CONFIG "SUNSHINE $attorb PREFILTER D D % ".
                     "is the sun shining on the satelite?\n";

        print CONFIG "ANG_DIST $attorb PREFILTER D D % ".
                     "degrees from nominal pointing\n";

        print CONFIG "RA $attorb PREFILTER D D % ".
                     "Right Ascension of satelite pointing\n";

        print CONFIG "DEC $attorb PREFILTER D D % ".
                     "Declination of satelite pointing\n";

        print CONFIG "ROLL $attorb PREFILTER D D % ".
                     "Roll angle of satelite pointing\n";

        print CONFIG "SAA $attorb PREFILTER D D % ".
                     "SAA derived from the orbit data\n";
                     
        print CONFIG "SAA_TIME $attorb PREFILTER D D % ".
                     "SAA time derived from the orbit data\n";
                     
        print CONFIG "COR_SAX $attorb PREFILTER D D % ".
                     "cutoff rigidity calculated as for the SAX\n";
                     
        print CONFIG "MCILWAIN_L $attorb PREFILTER D D % ".
                     "McIlwain L parameter for Earth's magnetic field\n";
                     
        print CONFIG "SUN_ANGLE $attorb PREFILTER D D % ".
                     "Angle between sun and pointing direction\n";

        print CONFIG "MOON_ANGLE $attorb PREFILTER D D % ".
                     "Angle between moon and pointing direction\n";

        print CONFIG "RAM_ANGLE $attorb PREFILTER D D % ".
                     "Angle between satelite motion and pointing\n";
    }

    #############################################################
    # collect information from the engineering HK
    #############################################################
    my $proto_config = $self->proctop()."/lists/engineering_hk_filters";
    open PROTO, "< $proto_config";
    
    my $last_inst="";
    my $hk;
    my @extnames;
    while(<PROTO>) {
        chomp;
        
        #####################
        # skip comments and blank lines
        #####################
        if(/^\s*#/ || /^\s*$/ ) { next }

        #################################
        # parse the current line
        #################################
        my ($inst, $mnemonic, $column, $apids, $comment ) = split /\s*\|\s*/, $_;

        if( $inst ne $last_inst) {
            ###############################
            # new instrument
            ###############################
            $last_inst = $inst;
            
            ###############################################
            # get the HK file for this instrument and
            # extract a list of the extensions it contains
            ################################################
            $hk = $filename->get('enhk', $inst, '', 0);
            $log->entry("Extracting filter quantites from $hk");

            @extnames=();
            if(-f $hk ) {
                @extnames = Util::FITSfile->new($hk)
                                        ->list_hdus();
            }
            
        } # end if this is a new instrument
        
        #######################################################
        # a mnemonic can be in multiple APIDS,
        # we need to determine which APIDS are present and
        # merge them if neccessary
        ######################################################
        $log->entry("Extracting $mnemonic, a.k.a. $column from APIDS $apids");
        my @apids = split /\s+/, $apids;
        my $list = Util::FITSlist->new();

        foreach my $apid (@apids) {
            my $ext = sprintf("hk%03xx001", $apid);

            if(grep {$_ eq $ext} @extnames ) {
                ##############################
                # we have data for this APID
                ##############################
                $list->add($hk, $ext);
            }
            
        }
        
        #################################################
        # if we don't have any data, then we will leave
        # this mnemonic out of the filter file.
        #################################################
        if($list->count() == 0 ) {
            $log->entry("No data for $mnemonic, a.k.a. $column");
            next;
        }
        
        #############################################
        # if we get here, then merge the data
        # don't worry, if we only have one APID, then
        # we just use that one extension
        #############################################
        my $tmp = $self->temp_file($mnemonic);
        my $merged=$list->merge($tmp, ' ', 'TIME', $mnemonic);

        ##########################################
        # sort the merged file if there was more
        # than one HK file
        ##########################################
        my $ext = $list->get_extension();
        if($merged eq $tmp) {
            ##########################
            # mulitple sources
            ##########################
            my $fits = Util::FITSfile->new($merged);
            $fits->sort("unique");
            $ext = 1; # merged data written to the first extension
        }

        ########################################################
        # write a row to the makefilter config file
        ########################################################
        print CONFIG "$mnemonic $merged $ext D D $column $comment\n";
    
    } # end of loop over mnemonics

    #######################################################
    # add XRT frame HK
    #######################################################
    my $framehead = $filename->get("hk", "xrt", "hd", 0);
    
    if ( -f $framehead ) {

        my %xrt_cols=(CCDTemp => "CCD Temperature",
                    Vod1 => "output drain voltage for left amp 1",
                    Vod2 => "output drain voltage for left amp 2",
                    Vrd1 => "reference voltage for amp 1",
                    Vrd2 => "reference voltage for amp 2",
                    Vsub => "substrate bias voltage",
                    Vbackjun => "back junction bias voltage",
                    Baselin1 => "Baseline voltage for Signal Chain A",
                    Baselin2 => "Baseline voltage for Signal Chain B",
		    ENDTIME => "",
		    PixGtULD => "");


        #####################################
        # loop over all the columns
        #####################################
        foreach my $column (keys %xrt_cols) {

            print CONFIG "$column $framehead FRAME D D % $xrt_cols{$column}\n";
        }
    } # end if we have a frame header file
    
    ###################################################
    # BAT DAP HK
    ###################################################
    my $dap = $filename->get("hk", "b", "dp", 0);
    if( -f $dap) {
        ################################################
        # here are all the columns we are interested in 
        ################################################
        my %cols=(BHVCtrl_Level  => "BAT High Voltage Control Levels",
                  BHVMon_Level   => "BAT High Voltage Monitor Levels",
                  DM_HK_Temp     => "BAT DM Side temperatures",
                  DM_HK_HVLkCurr => "BAT High Voltage Leakage Current");
                  
        ###########################################################
        # These are vector columns and we are just interested in 
        # the min and max values, so generate a CFITSIO extended
        # filename syntax statement which will calculate these.
        # At the same time we put entries into the makefilter config file
        # for these.
        # Some day we might add a mean or median
        ##########################################################
        my $minmax = $self->temp_file("dap_hk_min_max_values");
        my $spec="[col TIME;";
        foreach my $col (keys %cols) {
            foreach my $type ("MIN", "MAX") {
                my $newcol="${col}_${type}";

                $spec .= " $newcol(1I)=$type($col);";
                print CONFIG "$newcol $minmax DAP_HK D D % $cols{$col}\n";
            }
        }

        $spec .= "]";

        ######################################################
        # create a temporary file with the min/max values
        #####################################################
        my $fits = Util::FITSfile->new($dap, "DAP_HK");
        $fits->specs($spec);
    
        $log->entry("Extracting min/max values from $dap");
        $fits->copy($minmax);
        

    } # end if there is a BAT DAP HK file

    
    ###################################################
    # ACS bit flags
    ###################################################
    my $attitude = $filename->get("attitude", "s");
    
    if(-f $attitude) {
    
        ##################################################
        # Split the bit flags into separate columns in a
        # temporary file
        ##################################################
        my $flags1 = $self->temp_file("acs_flags1");
        my $flags2 = $self->temp_file("acs_flags2");

        my $acs_fits = Util::FITSfile->new($attitude, "ACS_DATA");

        $acs_fits->specs("[col TIME; ".
                              "TEN_ARCMIN(1B)=FLAGS .eq. b1xxxxxxx ; ".
                                 "SETTLED(1B)=FLAGS .eq. bx1xxxxxx]");

        $acs_fits->extract($flags1);
        
        $acs_fits->specs("[col TIME; ".
                                 "ACS_SAA(1B)=FLAGS .eq. bxx1xxxxx ; ".
                                "SAFEHOLD(1B)=FLAGS .eq. bxxx1xxxx]");

        $acs_fits->extract($flags2);
        
        ################################################################
        # now list the extracted columns to the config file
        ################################################################
        print CONFIG "TEN_ARCMIN $flags1 ACS_DATA D D % ACS reports within 10 arcmin of target\n";
        print CONFIG    "SETTLED $flags1 ACS_DATA D D % ACS reports settled on target\n";
        print CONFIG    "ACS_SAA $flags2 ACS_DATA D D % ACS reports in SAA\n";
        print CONFIG   "SAFEHOLD $flags2 ACS_DATA D D % ACS reports in SAFE\n";


    } # end if there is an attitude file



    #################################################
    # that's it for the config file, so close it up
    #################################################
    close CONFIG;

    ################################################
    # check if there is anything in the config file
    ################################################
    unless( -s $config ) {
        $log->entry("Nothing to put in the filter file");
        unlink $config;
        return;
    }


    ##############################
    # create the filter file
    ##############################
    my $makefilter = Util::HEAdas->new("makefilter");
    $log->entry("Creating filter file $filter");
    $makefilter->params({configure  => $config,
                         infileroot => "",
                         outfile    => $filter,
                         leapsec    => 'CALDB',
                         tprec      => 0.001,
                         mission    => "SWIFT",
                         clobber    => "yes"   });

    $makefilter->run();

    my $makefits = Util::FITSfile->new($filter);
    my @time = $makefits->cols('TIME')->table();
    my $first = Util::Date->new($time[0]);
    my $last  = Util::Date->new($time[-1]);

    $makefits->ext(0);
#    $makefits->keyword('TSTART', $time[0]);
#    $makefits->keyword('TSTOP', $time[-1]);
#    $makefits->keyword('TELAPSE', $time[-1]-$time[0]);
    $makefits->keyword('DATE-OBS', $first->date().'T'.$first->time() );
    $makefits->keyword('DATE-END', $last->date().'T'.$last->time() );

} # end of make_filter_file method


##############################################################################
# Make a number of plots of interesting quantities in the filter file
##############################################################################
sub make_plots {
    my $self = shift;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();

    ########################
    # get the filter file
    ########################
    my $filter = $filename->get("filter", "s");
    unless( -f $filter) { return; }

    my $fits_filter = Util::FITSfile->new($filter, 1);

    $log->entry("Plotting interesting things in the filter file");

    ######################################
    # set up for plotting
    ######################################
    my %yparm=(attd=>"RA DEC ROLL ANG_DIST",
               weel=>"FWHEEL"              ,
               bat1=>"BDIHKGRBFOLLOW BDIHKSLEWING BDIHKSAA BDIHKVALID",
               bat2=>"BDIHKNGOODDETS BCAHKDOINGBLK BBRHKTRIGNUM",
               saaa=>"ORBIT_SAA ACS_SAA BDIHKSAA"     );

    my $command_tmp=$self->temp_file("plot_commands");

    my $fplot=Util::Ftool->new("fplot")
                         ->params({infile=>$filter,
                                   xparm   => "TIME",
                                   rows    => "-",
                                   device  => "/cps",
                                   pltcmd => "\@$command_tmp",
                                   offset  => "yes"});

    my $plot_output="pgplot.ps";


    my $mode;
    foreach $mode (keys %yparm) {
        #################################################
        # make sure the filter file has all the columns
        # that we need for this plot
        #################################################
        my $good=1;
        foreach my $col (split /\s+/, $yparm{$mode}) {
            unless( $fits_filter->find_column($col)) {
                $log->entry("$filter is missing $col");
                $good=0;
                last;
            }
        }
        
        unless($good) {
            $log->entry("Skipping $mode plot");
            next;
        }

        #######################################
        # copy the plot command file and 
        # insert the file name for "MKF"
        #######################################
        my $template = $self->proctop()."/lists/filter_plot.$mode";
        open IN, "<$template";
        open OUT, ">$command_tmp";
        while (<IN>) {

            s/MKF/$filter/;
            print OUT;
        }
        close(OUT);
        close(IN);

        ##################################
        # make the plot
        ##################################
        unlink $plot_output;
        $fplot->params({yparm=>$yparm{$mode}})
              ->run();

        ############################################
        # save the plot with the correct file name
        ############################################
        if( -s $plot_output ) {
            rename $plot_output,
                   $filename->get("filterplot","s",$mode,0);
        } else {
            $log->error(2,"$plot_output not created");
        }
    } # end of loop over plot types


} # end of make_plots method


###########################################################################
# apply standard filters to the event data
###########################################################################
sub filter_data {

    my $self = shift;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();

	$self->xrtscreen;

	$self->uvotscreen;

} # end of filter_data'

sub filter_dataEAP {

    my $self = shift;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();
    
    #######################################
    # find the filter file
    #######################################
    my $filter = $filename->get("filter", "s");
    if( ! -f $filter ) {
        $log->entry("No filter file, so can't do screening");
        return;
    }


    ###########################
    # fire up the extractor
    ###########################
    my $extractor = Util::Extractor->new();
    $extractor->verbose(2);

    ########################################
    # now filter the data
    # this is going to need drastic changes
    ########################################
    my $inst;
    foreach $inst ("bat", "xrt", "uvot") {

        $log->entry("Filtering $inst files");

        ################################################
        # determine the modes which should be filtered
        ################################################
        my @modes=();
        if(   $inst eq "xrt" ) { @modes=("pcw?po","wtw?po", "puw?po", "lrw?po") }
        elsif($inst eq "uvot") { @modes=("*")    }
        elsif($inst eq "bat" ) { @modes=("evshps", "evshsl", "evshpo")    }

        ########################
        # loop over modes
        ########################
        my $mode;
        foreach $mode (@modes) {
            #################################################
            # determine what kid of files we should use as
            # input for filtering
            ################################################
            my $unfiltered_type="unfiltered";
            if($inst eq "xrt" && substr($mode,0,2) ne "pc" ) {
                ########################################################
                # for the XRT use the reconstructed event files unless
                # this is photon counting  mode
                ########################################################
                $unfiltered_type = "reconst";
            }

            #########################################
            # get the unfiltered files
            #########################################
            my @unfs = $filename->get($unfiltered_type, $inst, $mode, '*');
            
            #############################################
            # check if there are any files in this mode
            #############################################
            if(@unfs) {
                $log->entry("Mode $mode");
            } else {
                $log->entry("No unfiltered files for mode $mode\n");
                next;
            }

            ###############################################
            # get the mode name understood by the xselect
            # mission database
            ###############################################
            my $mdb_mode;
            if($inst eq "bat" ) { $mdb_mode = "EVENT" }
            elsif($inst eq "uvot" ) { $mdb_mode = "EVENT" }
            elsif($inst eq "xrt" ) {
                my $mo = substr($mode,0,2);
                if(   $mo eq "pc" ) { $mdb_mode="PHOTON"   }
                elsif($mo eq "wt" ) { $mdb_mode="WINDOWED" }
                elsif($mo eq "pu" ) { $mdb_mode="PILEDUP"  }
                elsif($mo eq "lr" ) { $mdb_mode="LOWRATE"  }
            
            }

            ####################################
            # set the mode and criteria to use
            ####################################
            $extractor->instrument($inst, $mdb_mode);

            ########################################
            # set the time filtering criteria
            ########################################
            if($inst eq "bat" ) {
                ##################################
                # BAT - no screening criteria yet
                # eventually just SAA probably
                ##################################
                $extractor->criteria("BDIHKNGOODDETS>1 && ".
                                     "BDIHKSAA==0 && ".
                                     "BHVMon_Level_MAX>99", $filter);

            } elsif($inst eq "xrt") {
                ######################################
                # XRT
                ######################################
                my $expression = $self->get_xrt_filter_expression();
                $extractor->criteria($expression, $filter);

            } elsif($inst eq "uvot") {
                $extractor->criteria("ELV>10 && ORBIT_SAA==0", $filter);
            }

            #######################################################
            # if there's no good time then there's no point doing
            # any filtering for this instrument
            ######################################################
            unless($extractor->has_good_time()) {
                $log->entry("No GTIs so skipping all $inst $mode files");
                next;
            }

            ###################################
            # region filter to use
            ###################################
            my $region="";
            if($inst eq "xrt" && $mode =~ /^pc/ ) {
                ##############################################
                # For XRT photon counting mode we screen out
                # the calibration sources
                ##############################################
                $region = $filename->fetch_cal("regstd", $inst, "");
            }

            $extractor->region($region);

            
            ################################################
            # get the row filter to use, if any
            # This is a CFITSIO extended filename syntax
            # row filtering statement which which will be
            # appended to the filename
            ################################################
            my $row_selection="";
            if($inst eq "xrt") {
                #######
                # XRT
                #######
                $row_selection = $self->get_xrt_row_selection($mdb_mode);
            } elsif($inst eq "uvot") {
                #############
                # UVOT
                #############
                $row_selection="QUALITY == 0";
            }

            $extractor->row_selection($row_selection);


            #############################
            # loop over unfiltered files
            #############################
            foreach my $unf (@unfs) {

                ##################################################
                # get the name of the corresponding filtered file
                ##################################################
                my $evt = $filename->corresponding($unfiltered_type, 'event', $unf);

                $log->entry("Extracting $evt from $unf");

                ###########################
                # run the extractor
                ###########################
                $extractor->infiles($unf)
                          ->outfile($evt, "event")
                          ->run();

                ##########################################
                # some special error checking
                ##########################################
                if(! $extractor->had_error() && ! -f $evt) {
                    $log->error(2, "Filtered file $evt not created");
                    next;
                }

                if($extractor->had_error() ) {
                    if( -f $evt) {
                        $log->entry("Deleting $evt");
                        unlink $evt;
                    }
                    next;
                }

                ####################################################
                # make sure the screened file has some events in it
                ####################################################
                my $nevents = Util::FITSfile->new($evt,"EVENTS")
                                            ->keyword("NAXIS2");
                if($nevents) {
                    $log->entry("$nevents filtered events");
                } else {
                    $log->entry("Deleteing $evt since it has $nevents events");
                    unlink $evt;
                    next;
                }

            } # end of loop over files
        } # end of loop over modes
    } # end of loop over instruments

} # end of filter_data method

############################################################################
#
############################################################################
sub read_xrt_range_file {
    my $self = shift;
    my $range_file=shift;
    my @extensions = @_;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();
    
        ##############################################
    # loop over the extensions of the range file
    ##############################################
    my @expression=();
    foreach my $ext (@extensions) {

        my $fits = Util::FITSfile->new($range_file, $ext);

        $fits->cols("PARNAME", "RANGE");
        my %range=$fits->table();

        ########################################
        # now loop over all the parameters
        ########################################
        foreach my $param (keys %range) {
            ##########################################
            # get the range string and
            # apply any aliases for the column names
            ##########################################
            my $range = $range{$param};


            #####################################
            # write the expression for this row
            #####################################
            if($range =~ /,/ ) {
                #############################
                # we have a min and or a max
                #############################
                my ($min, $max) = split /\s*,\s*/, $range;
                if($min ne "INDEF" ) { push @expression, ("$param >= $min") }
                if($max ne "INDEF" ) { push @expression, ("$param <= $max") }
            } else {
                ###################################################
                # no comma means force the param to a single value
                ###################################################
                push @expression, ("$param == $range");
            }
        }


    }

    return  join(" && ", @expression);

} # end of read_xrt_range_file method
    


############################################################################
# create a filter expression from an XRT-style parameter range file.
# This method remembers the expression it returned the first time it was
# run and will return the same value on subsequent calls without wasting
# any effort.
############################################################################
sub get_xrt_filter_expression {
    my $self = shift;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();

    if($self->{XRT_EXPRESSION} ) { return $self->{XRT_EXPRESSION} }

    ##################################
    # get the range file
    ##################################
    my $range_file = $filename->fetch_cal('rangefile', 'xrt');
    
    my $expression = $self->read_xrt_range_file($range_file, "ATTRANGE", "INSTRANGE");
    
    #####################################
    # alias "SAA" to "ORBIT_SAA"
    #####################################
    $expression =~ s/\sSAA\s/ ORBIT_SAA /;

    $self->{XRT_EXPRESSION} = $expression;

    return $self->{XRT_EXPRESSION};


} # end of make_expression_from_range_file method


############################################################################
#
############################################################################
sub get_xrt_row_selection {

    my $self = shift;
    my $datamode = shift;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();



    ################################################
    # find the event filter range calibration file
    ################################################
    my $range_file = $filename->fetch_cal('graderange', 'xrt');
    
    ###########################################
    # determine which extension we should read
    ###########################################
    my $ext="";
    if($datamode eq "LOWRATE" || $datamode eq "PILEDUP" ) { $ext = "PDRANGE" }
    elsif(                       $datamode eq "PHOTON"  ) { $ext = "PCRANGE" }
    elsif(                       $datamode eq "WINDOWED") { $ext = "WTRANGE" }
    else { $log->error(2, "Unknown DATAMODE $datamode")  }
    
    ##########################################
    # read the range file
    ##########################################
    my $expression = $self->read_xrt_range_file($range_file, $ext);

    return $expression;

} # end of get_xrt_row_filter method



sub uvotscreen {

    my $self = shift;

    my $log     = $self->log();
    my $filename= $self->filename();
    my $procpar = $self->procpar();
    my $jobpar  = $self->jobpar();

    my $attorb = $filename->get('attorb', 's');

    my $uscreen = Util::HEAdas->new('uvotscreen')
                     ->params({
                        attorbfile => $attorb,
                        aoexpr => 'ANG_DIST.lt.100',
                        evexpr => 'QUALITY==0',
                        badpixfile => 'CALDB'
                     })
           		->is_script(1);

    # iterate over UVOT unfiltered event files
    foreach my $unfFile ($filename->get('unfiltered', 'uvot', '*', '*')) {

        my $evtFile = $filename->corresponding('unfiltered', 'event', $unfFile);

        $log->entry("filtering $unfFile");

        $uscreen->params({
              infile => $unfFile,
              outfile => $evtFile,
              })
           ->run();
    }

} # end of uvotscreen


sub xrtscreen {

    my $self = shift;

    my $log     = $self->log();
    my $filename= $self->filename();
    my $procpar = $self->procpar();
    my $jobpar  = $self->jobpar();

    #######################################
    # find the filter file
    #######################################
    my $filter = $filename->get('filter', 's');
	if (not -f $filter) {
		$log->entry("No filter file, so can't do screening");
		return
	}

    ###################################
    # expression to use
    ###################################
    my $expression = $self->get_xrt_filter_expression();

	my $xrtscreen = Util::HEAdas->new('xrtscreen')
				->params({
					  createinstrgti => 'yes',
					  hkrangefile => 'CALDB',
					  evtrangefile => 'CALDB',
					  mkffile => $filter,
					  expr => $expression,
					 })
				->is_script(1);

    ################################################
    # determine the modes which should be filtered
    ################################################
    my @modes= qw(pcw?po wtw?po puw?po lrw?po);

	my $inst = 'xrt';
    ########################
    # loop over modes
    ########################
    foreach my $mode (@modes) {

		my $mo = substr($mode,0,2);

        #################################################
        # determine what kid of files we should use as
        # input for filtering
        ################################################
        my $unfiltered_type="unfiltered";
        if ($mo ne "pc") {
            ########################################################
            # for the XRT use the reconstructed event files unless
            # this is photon counting  mode
            ########################################################
            $unfiltered_type = "reconst";
        }

		{
			my $mdb_mode = '';

			if(   $mo eq "pc" ) { $mdb_mode="PHOTON"   }
			elsif($mo eq "wt" ) { $mdb_mode="WINDOWED" }
			elsif($mo eq "pu" ) { $mdb_mode="PILEDUP"  }
			elsif($mo eq "lr" ) { $mdb_mode="LOWRATE"  }

			my $gradexpr = $mdb_mode
					? $self->get_xrt_row_selection($mdb_mode)
					: 'none';

			# createattgti is only yes in PD unit test
			my $createattgti = $mo eq 'pd' ? 'yes' : 'no';

			# need to save a GTI file for each mode since it is required
			# for later xrtimage run
			my $gtifile = $filename->get('gti', 'xrt', $mo, 0);

			# evtscreen/gtiscreen is yes for all event modes, no for IM mode
			my $screen = $mo eq 'im' ? 'no' : 'yes';
			$xrtscreen->params({
					exprgrade => $gradexpr,
					createattgti => $createattgti,
					gtifile => $gtifile,
					gtiscreen => $screen,
					evtscreen => $screen,
					});
		}

        #########################################
        # get the unfiltered files
        #########################################
        my @unfs = $filename->get($unfiltered_type, $inst, $mode, '*');

        #############################################
        # check if there are any files in this mode
        #############################################
        if (@unfs) {
            $log->entry("Mode $mode");
        } else {
            $log->entry("No unfiltered files for mode $mode\n");
            next;
        }


        ###################################
        # region filter to use
        ###################################
        my $region="";
        if ($mo eq 'pc') {
            ##############################################
            # For XRT photon counting mode we screen out
            # the calibration sources
            ##############################################
            $region = $filename->fetch_cal("regstd", $inst, "");
        }


        #############################
        # loop over unfiltered files
        #############################
        foreach my $unf (@unfs) {

            ##################################################
            # get the name of the corresponding filtered file
            ##################################################
            my $evt = $filename->corresponding($unfiltered_type, 'event', $unf);

            $log->entry("Extracting $evt from $unf");

            ###########################
            # run xrtscreen
            ###########################
			$xrtscreen->params({
					infile => $unf,
					outfile => $evt,
					})
				->run;

			unlink(qw(xselect.hty defaults.def));

			if (not -f $evt) {
				$log->entry("Missing $evt");
				next;
			}

            ####################################################
            # make sure the screened file has some events in it
            ####################################################
            my $nevents = Util::FITSfile->new($evt,"EVENTS")
                                        ->keyword("NAXIS2");
            if($nevents) {
                $log->entry("$nevents filtered events");
            } else {
                $log->entry("Deleting $evt since it has $nevents events");
                unlink $evt;
                next;
            }

		} # end foreach unfiltered

	} # end foreach mode

} # end of xrtscreen



1;