Index: trunk/ippScripts/scripts/detrend_reject_exp.pl
===================================================================
--- trunk/ippScripts/scripts/detrend_reject_exp.pl	(revision 8510)
+++ trunk/ippScripts/scripts/detrend_reject_exp.pl	(revision 8510)
@@ -0,0 +1,150 @@
+#!/usr/bin/env perl
+
+use warnings;
+use strict;
+
+use IPC::Cmd qw( can_run run );
+use PS::IPP::Metadata::Config;
+use PS::IPP::Metadata::List qw( parse_md_list );
+use Statistics::Descriptive;
+use Data::Dumper;
+
+# Rejection threshold for the mean, in terms of the standard deviation
+# This measures how close it is to what's expected
+use constant REJECT_MEAN => {
+    'bias' => 0.01,		# Should be fairly flat; some CRs
+    'dark' => 0.05,		# Lots of CRs
+    'flat' => undef		# Can't define expected value (depends on exposure level)
+    };
+
+# Rejection threshold for the standard deviation, in terms of the standard deviation of the stdevs
+# This measures how much variation there is within the components of an exposure, compared to "typical"
+use constant REJECT_STDEV => {
+    'bias' => 3,		# Components should have the same degree of variation
+    'dark' => 3,		# Components should have the same degree of variation
+    'flat' => 3			# Components should have the same degree of variation
+    };
+
+# Rejection threshold for the mean of the stdev, in terms of the standard deviation of the mean of the stdevs
+# This measures how structured the images are, compared to "typical"
+use constant REJECT_MEAN_STDEV => {
+    'bias' => 3,		# All images should be equally structured
+    'dark' => 3,		# All images should be equally structured
+    'flat' => 3 		# All images should be equally structured
+    };
+
+# Parse command-line arguments
+die "Reject exposures based on the sample of exposures.\n\n" .
+    "Usage: $0 DET_ID ITER DETREND_TYPE\n\n"
+    if scalar @ARGV != 3;
+my $detId = shift @ARGV;	# Detrend ID
+my $iter = shift @ARGV;		# Iteration number
+my $detType = shift @ARGV;	# Detrend type
+
+# Look for programs we need
+my $missing_tools;
+my $dettool = can_run('dettool') or (warn "Can't find dettool" and $missing_tools = 1);
+die "Can't find required tools.\n" if $missing_tools;
+my $mdcParser = PS::IPP::Metadata::Config->new;	# Parser for metadata config files
+
+# Get list of component files
+my $exposures;			# Array of exposures
+{
+    my $command = "$dettool -residexp -det_id $detId -iter $iter"; # Command to run
+    my ( $success, $error_code, $full_buf, $stdout_buf, $stderr_buf ) =
+	run(command => $command, verbose => 1);
+    die "Unable to perform dettool -processed: $error_code\n" if not $success;
+    my $metadata = $mdcParser->parse(join "", @$stdout_buf); # Parsed metadata
+    $exposures = parse_md_list($metadata);
+}
+
+my @expIds;			# Array of exposure IDs
+my @means;			# Array of means
+my @stdevs;			# Array of stdevs
+my @meanStdevs;			# Array of mean stdevs
+foreach my $exposure (@$exposures) {
+    die "Unable to find exposure id.\n" if not defined $exposure->{exp_id};
+    die "Unable to find mean.\n" if not defined $exposure->{bg};
+    die "Unable to find stdev.\n" if not defined $exposure->{stdev};
+    die "Unable to find mean stdev.\n" if not defined $exposure->{bg_mean_stdev};
+    push @expIds, $exposure->{exp_id};
+    push @means, $exposure->{bg};
+    push @stdevs, $exposure->{bg_stdev};
+    push @meanStdevs, $exposure->{bg_mean_stdev};
+}
+my $meanStats = Statistics::Descriptive::Sparse->new(); # Statistics calculator
+$meanStats->add_data(@means);
+my $stdevStats = Statistics::Descriptive::Sparse->new(); # Statistics calculator
+$stdevStats->add_data(@stdevs);
+my $meanStdevStats = Statistics::Descriptive::Sparse->new(); # Statistics calculator
+$meanStdevStats->add_data(@meanStdevs);
+
+# Check the existence of the required rejection levels as a function of detrend type
+die "Unknown mean rejection level for detrend type $detType\n"
+    if not exists REJECT_MEAN->{$detType};
+die "Unknown stdev rejection level for detrend type $detType\n"
+    if not exists REJECT_STDEV->{$detType};
+die "Unknown mean stdev rejection level for detrend type $detType\n"
+    if not exists REJECT_MEAN_STDEV->{$detType};
+
+# Go through again to do rejection, and update the database for each exposure
+my $numRejected = 0;		# Number of exposures rejected
+for (my $i = 0; $i < scalar @means; $i++) {
+    my $expId = $expIds[$i];	# Exposure ID
+    my $command = "$dettool -updateresidexp -det_id $detId -iter $iter -exp_id $expId";	# Command to run
+    if ((defined REJECT_MEAN->{$detType} and
+	 ($means[$i] - $meanStats->mean()) / $meanStats->standard_deviation() > REJECT_MEAN->{$detType})
+	or
+	(defined REJECT_STDEV->{$detType} and
+	 ($stdevs[$i] - $stdevStats->mean()) / $stdevStats->standard_deviation() > REJECT_STDEV->{$detType})
+	or
+	(defined REJECT_MEAN_STDEV->{$detType} and
+	 ($meanStdevs[$i] - $meanStdevStats->mean()) / $meanStdevStats->standard_deviation >
+	 REJECT_MEAN_STDEV->{$detType})) {
+	$expId .= ' -reject';
+
+	### XXX: Need some way to know whether the exposure has already been rejected.
+	### rejection flag in $exposures?
+	$numRejected++;
+    }
+
+    # Update 
+    my ( $success, $error_code, $full_buf, $stdout_buf, $stderr_buf ) =
+	run(command => $command, verbose => 1);
+    die "Unable to perform dettool -updateresidexp: $error_code\n" if not $success;
+}
+
+# Decide if the current is sufficient to use as a master, and if we can stop iterating
+my $master = 1;			# This is good enough for a master
+my $stop = 1;			# Stop iterating
+# XXX: This probably isn't sufficient, but will do for now:
+if ($numRejected > 0) {
+    $master = 0;
+    $stop = 0;
+}
+
+# Put the result into the database
+{
+    my $command = "$dettool -adddetrunsummary -det_id $detId -iter $iter " .
+	"-bg " . $meanStats->mean() . " -bg_stdev " . $stdevStats->mean() .
+	" -bg_mean_stdev " . $meanStdevStats->mean();
+    $command .= " -reject" if not $master;
+    my ( $success, $error_code, $full_buf, $stdout_buf, $stderr_buf ) =
+	run(command => $command, verbose => 1);
+    die "Unable to perform dettool -adddetrunsummary: $error_code\n" if not $success;
+}
+
+# Re-run processing if required
+if (not $stop) {
+    my $command = "$dettool -updatedetrun -det_id $detId -iter $iter";
+    if ($stop) {
+	$command .= ' -stop';
+    } else {
+	$command .= ' -go';
+    }
+    my ( $success, $error_code, $full_buf, $stdout_buf, $stderr_buf ) =
+	run(command => $command, verbose => 1);
+    die "Unable to perform dettool -updatedetrun: $error_code\n" if not $success;
+}
+
+__END__