Index: trunk/ippScripts/scripts/detrend_reject_exp.pl
===================================================================
--- trunk/ippScripts/scripts/detrend_reject_exp.pl	(revision 13962)
+++ trunk/ippScripts/scripts/detrend_reject_exp.pl	(revision 14009)
@@ -47,12 +47,33 @@
 
 pod2usage( -msg => "Unknown option: @ARGV", -exitval => 2 ) if @ARGV;
-pod2usage(
-    -msg => "Required options: --det_id --iteration --det_type --camera",
-    -exitval => 3,
-) unless defined $det_id
+pod2usage( -msg => "Required options: --det_id --iteration --det_type --camera",
+	   -exitval => 3)
+    unless defined $det_id
     and defined $iter
     and defined $det_type
     and defined $camera;
 
+# values to extract from output metadata and the stats to calculate
+# XXX -bg_mean_stdev should take rms of bg_mean_stdev if bg_mean_stdev != 0 (A)
+# XXX -bg_mean_stdev should take stdev of bg_mean if bg_mean_stdev == 0     (B)
+# XXX  (A) if imfile.Ncomp > 1, (B) if imfile.Ncomp == 1
+my $STATS = 
+   [   
+       #          KEYWORD                 STATISTIC          CHIPTOOL FLAG
+       { name => "bg",             type => "mean",  flag => "-bg" },
+       { name => "bg_mean_stdev",  type => "stdev", flag => "-bg_mean_stdev" },
+       { name => "bg_stdev",       type => "rms",   flag => "-bg_stdev" },
+   ];
+
+my $REJSTATS = 
+   [   
+       #          KEYWORD                 STATISTIC          CHIPTOOL FLAG
+       { name => "bg",             type => "mean",   flag => "ensMeanMean" },
+       { name => "bg",             type => "stdev",  flag => "ensMeanStdev" },
+       { name => "bg_mean_stdev",  type => "mean",   flag => "ensMeanStdevMean" },
+       { name => "bg_mean_stdev",  type => "stdev",  flag => "ensMeanStdevStdev" },
+       { name => "bg_stdev",       type => "mean",   flag => "ensStdevMean" },
+       { name => "bg_stdev",       type => "stdev",  flag => "ensStdevStdev" },
+   ];
 
 # Look for programs we need
@@ -64,10 +85,11 @@
 }
 
-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 $det_id -iteration $iter"; # Command to run
+    # dettool command to select exp data for this det_run
+    my $command = "$dettool -residexp";
+    $command .= " -det_id $det_id";
+    $command .= " -iteration $iter";
     $command .= " -dbname $dbname" if defined $dbname;
     my ( $success, $error_code, $full_buf, $stdout_buf, $stderr_buf ) =
@@ -78,46 +100,42 @@
     }
 
+    # Parse the stdout buffer into a metadata
+    my $mdcParser = PS::IPP::Metadata::Config->new;	# Parser for metadata config files
     my $metadata = $mdcParser->parse(join "", @$stdout_buf) or
 	&my_die("Unable to parse metadata config doc", $det_id, $iter, $PS_EXIT_PROG_ERROR);
-    # XXX need to test that the result list is 0 length and return with a different error
+
+    # parse the file info in the metadata
     $exposures = parse_md_list($metadata) or
 	&my_die("Unable to parse metadata list", $det_id, $iter, $PS_EXIT_PROG_ERROR);
+
+    # Parse the statistics on the residual image
+    my $stats = PS::IPP::Metadata::Stats->new($STATS); # Stats parser
+    $stats->parse($metadata) or &my_die("Unable to find all values in statistics output.", $det_id, $iter, $PS_EXIT_PROG_ERROR);
+
+    # Parse the statistics for rejections
+    my $rejstats = PS::IPP::Metadata::Stats->new($REJSTATS); # Stats parser
+    $rejstats->parse($metadata) or &my_die("Unable to find all values in statistics output.", $det_id, $iter, $PS_EXIT_PROG_ERROR);
 }
 
 my @expTags;			# Array of exposure IDs
-my @means;			# Array of means
-my @stdevs;			# Array of standard deviations
-my @meanStdevs;			# Array of standard deviations of the mean (normalised by the mean)
-my @variances;			# Array of variances
 my @accept;			# Array of accept flags
 my @include;			# Array of include flags
 foreach my $exposure (@$exposures) {
     &my_die("Unable to find exposure id.\n", $det_id, $iter, $PS_EXIT_SYS_ERROR) unless defined $exposure->{exp_tag};
-    &my_die("Unable to find mean.\n", $det_id, $iter, $PS_EXIT_SYS_ERROR) unless defined $exposure->{bg};
-    &my_die("Unable to find stdev.\n", $det_id, $iter, $PS_EXIT_SYS_ERROR) unless defined $exposure->{bg_stdev};
-    &my_die("Unable to find mean stdev.\n", $det_id, $iter, $PS_EXIT_SYS_ERROR) unless defined $exposure->{bg_mean_stdev};
     &my_die("Unable to find accept.\n", $det_id, $iter, $PS_EXIT_SYS_ERROR) unless defined $exposure->{accept};
     &my_die("Unable to find include.\n", $det_id, $iter, $PS_EXIT_SYS_ERROR) unless defined $exposure->{include};
-    push @expTags, $exposure->{exp_tag};
-    push @means, $exposure->{bg};
-    push @stdevs, $exposure->{bg_stdev};
-    # XXX why are we using mean_stdev / bd?
-    # push @meanStdevs, $exposure->{bg_mean_stdev} / $exposure->{bg};
-    push @meanStdevs, $exposure->{bg_mean_stdev};
-    push @variances, $exposure->{bg_stdev}**2;
-    push @accept, $exposure->{accept};
-    push @include, $exposure->{include};
-}
-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);
+}
+
+my $ensMeanMean       = &STATS_value_for_flag ($REJSTATS, "ensMeanMean");      
+my $ensMeanStdev      = &STATS_value_for_flag ($REJSTATS, "ensMeanStdev");     
+my $ensMeanStdevMean  = &STATS_value_for_flag ($REJSTATS, "ensMeanStdevMean"); 
+my $ensMeanStdevStdev = &STATS_value_for_flag ($REJSTATS, "ensMeanStdevStdev");
+my $ensStdevMean      = &STATS_value_for_flag ($REJSTATS, "ensStdevMean");     
+my $ensStdevStdev     = &STATS_value_for_flag ($REJSTATS, "ensStdevStdev");    
 
 $ipprc->define_camera($camera);
 # Rejection thresholds
-my $reject_mean = rejection_limit( 'ENSEMBLE.MEAN', $det_type, $filter );
-my $reject_stdev = rejection_limit( 'ENSEMBLE.STDEV', $det_type, $filter );
+my $reject_mean      = rejection_limit( 'ENSEMBLE.MEAN',      $det_type, $filter );
+my $reject_stdev     = rejection_limit( 'ENSEMBLE.STDEV',     $det_type, $filter );
 my $reject_meanstdev = rejection_limit( 'ENSEMBLE.MEANSTDEV', $det_type, $filter );
 
@@ -129,6 +147,5 @@
 unless ($no_op) {
     $logFile = $ipprc->file_create_append( $logName );
-    print $logFile "Ensemble mean " . $meanStats->mean() . " +/- " . $meanStats->standard_deviation .
-	", stdev " . $stdevStats->mean() . " +/- " . $stdevStats->standard_deviation() . "\n\n";
+    print $logFile "Ensemble mean $meanEnsemble +/- $stdevEnsemble, stdev $meanStdevEnsemble\n\n";
 }
 
@@ -136,15 +153,24 @@
 my $numChanges = 0;		# Number of exposures with changed status
 my $numReject = 0;		# Number of exposures rejected
-my $command;
-my $reject;
-my $expTag;
-for (my $i = 0; $i < scalar @means; $i++) {
-    $expTag = $expTags[$i];	# Exposure ID
-    $command = "$dettool -updateresidexp -det_id $det_id -iteration $iter -exp_tag $expTag"; # Command to run
+
+for (my $i = 0; $i < scalar @$exposures; $i++) {
+    my $file      = $files[$i];
+    my $mean      = $file->{bg};	# Mean for this exposure
+    my $stdev     = $file->{bg_stdev}; # Stdev for this exposure
+    my $meanStdev = $file->{bg_mean_stdev}; # Stdev of Means for this exposure
+
+    my $expTag  = $exposure->{exp_tag};
+    my $accept  = $exposure->{accept};
+    my $include = $exposure->{include};
+
+    my $reject = 0;		# Reject this exposure?
+
+    my $command = "$dettool -updateresidexp";
+    $command .= " -det_id $det_id";
+    $command .= " -iteration $iter";
+    $command .= " -exp_tag $expTag";
     $command .= " -dbname $dbname" if defined $dbname;
     
-    $reject = 0;		# Reject this exposure?
-
-    if (not $accept[$i]) {
+    if (not $accept) {
 	# Rejected this at an earlier stage
 	unless ($no_op) {
@@ -165,7 +191,7 @@
     }
     
-    if ($reject_mean > 0 and defined $meanStats->standard_deviation() ) {
-	my $dMean = abs($means[$i] - $meanStats->mean()) ;
-	if ($dMean > ($reject_mean * $meanStats->standard_deviation())) {
+    if ($reject_mean > 0 and $ensMeanStdev > 0) {
+	my $delta = abs($mean - $ensMeanMean);
+	if ($delta > ($reject_mean * $ensMeanStdev)) {
 	    print $logFile "Rejecting $expTag based on ensemble mean value: ";
 	    $reject = 1;
@@ -174,12 +200,12 @@
 	    print $logFile "$expTag OK against ensemble mean: ";
 	}
-	print $logFile "$means[$i] --> $dMean vs " . $reject_mean * $meanStats->standard_deviation() . "\n";
+	print $logFile "$mean --> $delta vs " . $reject_mean * $ensMeanStdev . "\n";
     } else {
 	print $logFile "No rejection of $expTag for ensemble mean\n";
     }
 
-    if ($reject_stdev > 0 and defined $stdevStats->standard_deviation() > 0) {
-	my $dMean = abs($stdevs[$i] - $stdevStats->mean());
-	if ($dMean > ($reject_stdev * $stdevStats->standard_deviation())) {
+    if ($reject_stdev > 0 and $ensStdevStdev > 0) {
+	my $delta = abs($stdev - $ensStdevMean);
+	if ($delta > ($reject_stdev * $ensStdevStdev)) {
 	    print $logFile "Rejecting $expTag based on ensemble stdev: ";
 	    $reject = 1;
@@ -188,12 +214,12 @@
 	    print $logFile "$expTag OK against ensemble stdev: ";
 	}
-	print $logFile "$stdevs[$i] --> $dMean sigma vs " . $reject_stdev  * $stdevStats->standard_deviation() . "\n";
+	print $logFile "$stdev --> $delta sigma vs " . $reject_stdev * $ensStdevStdev . "\n";
     } else {
 	print $logFile "No rejection of $expTag for ensemble stdev\n";
     }
     
-    if ($reject_meanstdev > 0 and defined $meanStdevStats->standard_deviation() > 0) {
-	my $dMean = abs($meanStdevs[$i] - $meanStdevStats->mean());
-	if ($dMean > ($reject_meanstdev * $meanStdevStats->standard_deviation())) {
+    if ($reject_meanstdev > 0 and $ensMeanStdevStdev > 0) {
+	my $delta = abs($meanStdev - $ensMeanStdevMean);
+	if ($delta > ($reject_meanstdev * $ensMeanStdevStdev)) {
 	    print $logFile "Rejecting $expTag based on ensemble mean stdev: ";
 	    $reject = 1;
@@ -202,5 +228,5 @@
 	    print $logFile "$expTag OK against ensemble mean stdev: ";
 	}
-	print $logFile "$meanStdevs[$i] --> $dMean sigma vs " . $reject_meanstdev . "\n";
+	print $logFile "$meanStdev --> $delta sigma vs " . $reject_meanstdev * $ensMeanStdevStdev. "\n";
     } else {
 	print $logFile "No rejection of $expTag for ensemble mean stdev\n";
@@ -214,5 +240,5 @@
     
     # Check for status changes
-    if ((not $include[$i] and not $reject) or ($include[$i] and $reject)) {
+    if ((not $include and not $reject) or ($include and $reject)) {
 	unless ($no_op) {
 	    print $logFile "Status of $expTag has changed.\n";
@@ -253,49 +279,44 @@
 }
 
+my $command = "$dettool -adddetrunsummary";
+$command .= " -det_id $det_id";
+$command .= " -iteration $iter";
+$command .= " -accept" if $master;
+$command .= " -dbname $dbname" if defined $dbname;
+
+# add in the elements from the selected stats above
+foreach my $entry (@$STATS) {
+    my $value = $entry->{value};
+    my $flag = $entry->{flag};
+    $command .= " $flag $value";
+}
+
 # Put results into the database
 unless ($no_update) {
-
-    {
-	# Add the summary
-	my $varianceStats = Statistics::Descriptive::Sparse->new(); # Statistics calculator
-	$varianceStats->add_data(@variances);
-	
-	my $bg = ($meanStats->mean() or 'NAN');
-	my $bg_stdev = (sqrt( $varianceStats->mean() ) or 'NAN');
-	my $bg_mean_stdev = ($meanStats->standard_deviation() or 'NAN');
-	
-	my $command = "$dettool -adddetrunsummary -det_id $det_id -iteration $iter";
-	$command .= " -bg $bg -bg_stdev $bg_stdev -bg_mean_stdev $bg_mean_stdev";
-	$command .= " -accept" if $master;
-	$command .= " -dbname $dbname" if defined $dbname;
-	
-	my ( $success, $error_code, $full_buf, $stdout_buf, $stderr_buf ) =
-	    run(command => $command, verbose => 1);
-	unless ($success) {
-	    $error_code = (($error_code >> 8) or $PS_EXIT_PROG_ERROR);
-	    warn("Unable to perform dettool -adddetrunsummary: $error_code");
-	    exit($error_code);
-	}
-    }
-
-    # Re-run processing if required
-    {
-	my $command = "$dettool -updatedetrun -det_id $det_id";
-	if ($stop) {
-	    $command .= ' -state stop';
-	} else {
-	    $command .= ' -again';
-	}
-	
-	$command .= " -dbname $dbname" if defined $dbname;
-	
-	my ( $success, $error_code, $full_buf, $stdout_buf, $stderr_buf ) =
-	    run(command => $command, verbose => 1);
-	unless ($success) {
-	    $error_code = (($error_code >> 8) or $PS_EXIT_PROG_ERROR);
-	    warn("Unable to perform dettool -updatedetrun: $error_code");
-	    exit($error_code);
-	}
-    }
+    my ( $success, $error_code, $full_buf, $stdout_buf, $stderr_buf ) =
+	run(command => $command, verbose => 1);
+    unless ($success) {
+	$error_code = (($error_code >> 8) or $PS_EXIT_PROG_ERROR);
+	warn("Unable to perform dettool -adddetrunsummary: $error_code");
+	exit($error_code);
+    }
+
+    $command = "$dettool -updatedetrun -det_id $det_id";
+    if ($stop) {
+	$command .= ' -state stop';
+    } else {
+	$command .= ' -again';
+    }
+    $command .= " -dbname $dbname" if defined $dbname;
+
+    my ( $success, $error_code, $full_buf, $stdout_buf, $stderr_buf ) =
+	run(command => $command, verbose => 1);
+    unless ($success) {
+	$error_code = (($error_code >> 8) or $PS_EXIT_PROG_ERROR);
+	warn("Unable to perform dettool -updatedetrun: $error_code");
+	exit($error_code);
+    }
+} else {
+    print "skipping command: $command\n";
 }
 
@@ -317,5 +338,8 @@
     carp($msg);
     if ($det_id and $iter and not $no_update) {
-	my $command = "$dettool -adddetrunsummary -det_id $det_id -iteration $iter -code $exit_code";
+	my $command = "$dettool -adddetrunsummary";
+	$command .= " -det_id $det_id";
+	$command .= " -iteration $iter";
+	$command .= " -code $exit_code";
 	$command .= " -dbname $dbname" if defined $dbname;
 ###        system ($command);
@@ -323,6 +347,4 @@
     exit $exit_code;
 }
-
-
 
 # Retrieve the requested rejection limit, dying unless extant
@@ -342,4 +364,16 @@
 }
 
+sub STATS_value_for_flag 
+{
+    my $STATS = shift;
+    my $flag  = shift;
+
+    foreach my $entry (@$STATS) {
+	if ($flag eq $entry->{flag}) {
+	    return $entry->{value};
+	}
+    }
+    return 'NAN';
+}
 
 __END__
