Changeset 17943 for trunk/ippScripts/scripts/detrend_reject_imfile.pl
- Timestamp:
- Jun 5, 2008, 12:09:48 PM (18 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/ippScripts/scripts/detrend_reject_imfile.pl
r17671 r17943 3 3 # this program has two jobs: 4 4 # 1: for the given exp_tag, generate the binned & mosaiced JPEG images 5 # 2: examine the collection of per-imfile statistics and reject. At the moment, 5 # 2: examine the collection of per-imfile statistics and reject. At the moment, 6 6 # this program (and the database) only allows rejection at the exposure level, 7 # not at the imfile level. 7 # not at the imfile level. 8 8 9 9 use Carp; … … 32 32 33 33 use Getopt::Long qw( GetOptions :config auto_help auto_version gnu_getopt ); # option parsing 34 use Pod::Usage qw( pod2usage ); 34 use Pod::Usage qw( pod2usage ); 35 35 36 36 # parse the command-line options … … 38 38 $verbose, $no_update, $no_op, $save_temps ); 39 39 GetOptions( 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 40 'det_id|d=s' => \$det_id, 41 'iteration=s' => \$iter, 42 'exp_id|e=s' => \$exp_id, 43 'exp_tag|=s' => \$exp_tag, 44 'det_type|t=s' => \$det_type, 45 'camera=s' => \$camera, 46 'outroot|w=s' => \$outroot, # output file base name 47 'filter=s' => \$filter, 48 'reject' => \$reject, 49 'dbname|d=s' => \$dbname, # Database name 50 'reduction|=s' => \$reduction, 51 'verbose' => \$verbose, # Print to stdout 52 'no-update' => \$no_update, 53 'no-op' => \$no_op, 54 'save-temps' => \$save_temps, # Save temporary files? 55 ) or pod2usage( 2 ); 56 56 57 57 pod2usage( -msg => "Unknown option: @ARGV", -exitval => 2 ) if @ARGV; 58 58 pod2usage( -msg => "Required options: --det_id --iteration --exp_id --exp_tag --det_type --camera --outroot", 59 -exitval => 3) unless 60 defined $det_id and 61 defined $iter and 62 defined $exp_id and 63 defined $exp_tag and 64 defined $det_type and 59 -exitval => 3) unless 60 defined $det_id and 61 defined $iter and 62 defined $exp_id and 63 defined $exp_tag and 64 defined $det_type and 65 65 defined $camera and 66 66 defined $outroot; … … 82 82 # XXX -bg_mean_stdev should take stdev of bg_mean if bg_mean_stdev == 0 (B) 83 83 # XXX (A) if imfile.Ncomp > 1, (B) if imfile.Ncomp == 1 84 my $STATS = 84 my $STATS = 85 85 [ 86 87 88 89 { name => "bg_stdev", type => "rms", flag => "-bg_stdev", dtype => "float" }, 90 { name => "bg_skewness", type => "mean", flag => "-bg_skewness", dtype => "float" }, 91 { name => "bg_kurtosis", type => "mean", flag => "-bg_kurtosis", dtype => "float" }, 92 93 94 95 96 { name => "user_1", type => "mean", flag => "-user_1",dtype => "float" }, # fringe residual97 { name => "user_2", type => "rms", flag => "-user_2",dtype => "float" }, # fringe residual98 { name => "user_1", type => "stdev", flag => "-user_3",dtype => "float" }, # fringe residual99 86 # KEYWORD STATISTIC CHIPTOOL FLAG 87 { name => "bg", type => "mean", flag => "-bg", dtype => "float" }, 88 { name => "bg_mean_stdev", type => "stdev", flag => "-bg_mean_stdev", dtype => "float" }, 89 { name => "bg_stdev", type => "rms", flag => "-bg_stdev", dtype => "float" }, 90 { name => "bg_skewness", type => "mean", flag => "-bg_skewness", dtype => "float" }, 91 { name => "bg_kurtosis", type => "mean", flag => "-bg_kurtosis", dtype => "float" }, 92 { name => "bin_stdev", type => "rms", flag => "-bin_stdev", dtype => "float" }, 93 { name => "fringe_0", type => "mean", flag => "-fringe_0", dtype => "float" }, 94 { name => "fringe_1", type => "rms", flag => "-fringe_1", dtype => "float" }, 95 { name => "fringe_0", type => "stdev", flag => "-fringe_2", dtype => "float" }, 96 { name => "user_1", type => "mean", flag => "-user_1", dtype => "float" }, # fringe residual 97 { name => "user_2", type => "rms", flag => "-user_2", dtype => "float" }, # fringe residual 98 { name => "user_1", type => "stdev", flag => "-user_3", dtype => "float" }, # fringe residual 99 ]; 100 100 my $stats = PS::IPP::Metadata::Stats->new($STATS); # Stats parser 101 101 … … 104 104 my $dettool = can_run('dettool') or (warn "Can't find dettool" and $missing_tools = 1); 105 105 my $ppImage = can_run('ppImage') or (warn "Can't find ppImage" and $missing_tools = 1); 106 if ($missing_tools) { 106 if ($missing_tools) { 107 107 warn("Can't find required tools."); 108 exit($PS_EXIT_CONFIG_ERROR); 108 exit($PS_EXIT_CONFIG_ERROR); 109 109 } 110 110 … … 119 119 $command .= " -dbname $dbname" if defined $dbname; 120 120 ( $success, $error_code, $full_buf, $stdout_buf, $stderr_buf ) = 121 121 run(command => $command, verbose => $verbose); 122 122 unless ($success) { 123 124 125 123 $error_code = (($error_code >> 8) or $PS_EXIT_PROG_ERROR); 124 warn("Unable to perform dettool -residimfile: $error_code\n"); 125 exit($error_code); 126 126 } 127 127 # XXX report an error message if stdout_buf is empty 128 128 129 129 # Parse the stdout buffer into a metadata 130 my $mdcParser = PS::IPP::Metadata::Config->new; 130 my $mdcParser = PS::IPP::Metadata::Config->new; 131 131 my $metadata = $mdcParser->parse(join "", @$stdout_buf) or 132 132 &my_die("Unable to parse metadata config doc", $det_id, $iter, $exp_id, $PS_EXIT_PROG_ERROR); 133 133 134 134 # parse the file info in the metadata 135 135 $files = parse_md_list($metadata) or 136 136 &my_die("Unable to parse metadata list", $det_id, $iter, $exp_id, $PS_EXIT_PROG_ERROR); 137 137 138 138 # Parse the statistics on the residual image … … 172 172 # Make the jpeg for binning 1 173 173 $command = "$ppImage -list $list1Name $outroot -recipe PPIMAGE $recipe1"; # Command to run 174 $command .= " -dbname $dbname" if defined $dbname; 174 175 ( $success, $error_code, $full_buf, $stdout_buf, $stderr_buf ) = 175 176 run(command => $command, verbose => $verbose); 176 177 unless ($success) { 177 178 178 $error_code = (($error_code >> 8) or $PS_EXIT_PROG_ERROR); 179 &my_die("Unable to run ppImage: $error_code", $det_id, $iter, $exp_id, $error_code); 179 180 } 180 181 &my_die("Unable to find expected output file: $jpeg1Name", $det_id, $iter, $exp_id, $PS_EXIT_SYS_ERROR) unless -f $ipprc->file_resolve($jpeg1Name); 181 182 182 183 # Make the jpeg for binning 2 183 184 $command = "$ppImage -list $list2Name $outroot -recipe PPIMAGE $recipe2"; # Command to run 185 $command .= " -dbname $dbname" if defined $dbname; 184 186 ( $success, $error_code, $full_buf, $stdout_buf, $stderr_buf ) = 185 187 run(command => $command, verbose => $verbose); 186 188 unless ($success) { 187 188 189 $error_code = (($error_code >> 8) or $PS_EXIT_PROG_ERROR); 190 &my_die("Unable to run ppImage: $error_code", $det_id, $iter, $exp_id, $error_code); 189 191 } 190 192 &my_die("Unable to find expected output file: $jpeg2Name", $det_id, $iter, $exp_id, $PS_EXIT_SYS_ERROR) unless -f $ipprc->file_resolve($jpeg2Name); … … 217 219 foreach my $file (@$files) { 218 220 my $name = $file->{class_id}; 219 my $mean = $file->{bg}; 221 my $mean = $file->{bg}; # Mean for this imfile 220 222 my $stdev = $file->{bg_stdev}; # Stdev for this imfile 221 223 my $skewness = $file->{bg_skewness}; # Skewness for this imfile … … 227 229 my $flux; 228 230 if ($file->{exp_time} == 0.0) { 229 230 } else { 231 231 $flux = $mean; 232 } else { 233 $flux = $mean / $file->{exp_time}; 232 234 } 233 235 push @fluxes, $flux; … … 243 245 # 0.0 for a bias). in general, use with ADDITIVE components 244 246 if ($reject_imfile_mean > 0) { 245 246 247 248 249 250 251 247 if (abs($mean) > $reject_imfile_mean * $stdev) { 248 print $logFile "Rejecting exposure based on bad imfile mean for $name: "; 249 $reject = 1; 250 } else { 251 print $logFile "Imfile mean for $name meets requirements: "; 252 } 253 print $logFile "$mean vs $reject_imfile_mean * $stdev\n"; 252 254 } else { 253 255 print $logFile "No rejection on imfile mean for $name\n"; 254 256 } 255 257 … … 260 262 # for a bias). in general, use with ADDITIVE components 261 263 if ($reject_imfile_flux > 0) { 262 263 264 265 266 267 268 264 if (abs($flux) > $reject_imfile_flux) { 265 print $logFile "Rejecting exposure based on bad imfile flux for $name: "; 266 $reject = 1; 267 } else { 268 print $logFile "Imfile flux for $name meets requirements: "; 269 } 270 print $logFile "$flux vs $reject_imfile_flux\n"; 269 271 } else { 270 272 print $logFile "No rejection on imfile flux for $name\n"; 271 273 } 272 274 … … 274 276 # larger than the limit 275 277 if ($reject_imfile_stdev > 0) { 276 277 278 279 280 281 282 283 284 } else { 285 278 if ($stdev > $reject_imfile_stdev) { 279 print $logFile "Rejecting exposure based on bad imfile stdev for $name: "; 280 $reject = 1; 281 } else { 282 print $logFile "Imfile stdev for $name meets requirements: "; 283 } 284 print $logFile "$stdev vs $reject_imfile_stdev\n"; 285 286 } else { 287 print $logFile "No rejection on imfile stdev for $name\n"; 286 288 } 287 289 … … 289 291 # larger than the limit 290 292 if ($reject_imfile_skewness > 0) { 291 292 293 294 295 296 297 298 299 } else { 300 293 if ($stdev > $reject_imfile_skewness) { 294 print $logFile "Rejecting exposure based on bad imfile skewness for $name: "; 295 $reject = 1; 296 } else { 297 print $logFile "Imfile skewness for $name meets requirements: "; 298 } 299 print $logFile "$skewness vs $reject_imfile_skewness\n"; 300 301 } else { 302 print $logFile "No rejection on imfile skewness for $name\n"; 301 303 } 302 304 … … 304 306 # larger than the limit 305 307 if ($reject_imfile_kurtosis > 0) { 306 307 308 309 310 311 312 313 314 } else { 315 308 if ($stdev > $reject_imfile_kurtosis) { 309 print $logFile "Rejecting exposure based on bad imfile kurtosis for $name: "; 310 $reject = 1; 311 } else { 312 print $logFile "Imfile kurtosis for $name meets requirements: "; 313 } 314 print $logFile "$kurtosis vs $reject_imfile_kurtosis\n"; 315 316 } else { 317 print $logFile "No rejection on imfile kurtosis for $name\n"; 316 318 } 317 319 … … 319 321 # component means is larger than the limit 320 322 if ($reject_imfile_meanstdev > 0) { 321 322 323 324 325 326 327 328 } else { 329 323 if ($meanStdev > $reject_imfile_meanstdev) { 324 print $logFile "Rejecting exposure based on bad imfile mean stdev for $name: "; 325 $reject = 1; 326 } else { 327 print $logFile "Imfile mean stdev for $name meets requirements: "; 328 } 329 print $logFile "$meanStdev vs $reject_imfile_meanstdev\n"; 330 } else { 331 print $logFile "No rejection on imfile mean stdev for $name\n"; 330 332 } 331 333 … … 333 335 # the mean counts / total pixel variance) of the imfile component 334 336 # means are less than the limit. this test is sensible for images 335 # which have finite residual flux such as a flat-field image. 337 # which have finite residual flux such as a flat-field image. 336 338 if ($reject_imfile_snr > 0) { 337 338 339 340 341 342 343 344 } else { 345 339 if ($mean < $stdev * $reject_imfile_snr) { 340 print $logFile "Rejecting exposure based on bad imfile S/N for $name: "; 341 $reject = 1; 342 } else { 343 print $logFile "Imfile S/N for $name meets requirements: "; 344 } 345 print $logFile "mean: $mean vs stdev*SNlimit: " . $stdev * $reject_imfile_snr . "\n"; 346 } else { 347 print $logFile "No rejection on imfile S/N for $name\n"; 346 348 } 347 349 … … 351 353 # have finite residual flux such as a flat-field image. 352 354 if ($reject_imfile_bin_stdev > 0) { 353 354 355 356 357 358 359 360 } else { 361 355 if ($binStdev > $reject_imfile_bin_stdev) { 356 print $logFile "Rejecting exposure based on bad imfile binned stdev for $name: "; 357 $reject = 1; 358 } else { 359 print $logFile "Imfile binned stdev for $name meets requirements: "; 360 } 361 print $logFile "$binStdev vs $reject_imfile_bin_stdev\n"; 362 } else { 363 print $logFile "No rejection on imfile binned stdev for $name\n"; 362 364 } 363 365 if ($reject_imfile_bin_snr > 0) { 364 365 366 367 368 369 370 371 } else { 372 366 if ($mean < $binStdev * $reject_imfile_bin_snr) { 367 print $logFile "Rejecting exposure based on bad imfile binned S/N for $name: "; 368 $reject = 1; 369 } else { 370 print $logFile "Imfile binned S/N for $name meets requirements: "; 371 } 372 print $logFile "mean: $mean vs binStdev*SNlimit: " . $binStdev * $reject_imfile_bin_snr . "\n"; 373 } else { 374 print $logFile "No rejection on imfile binned S/N for $name\n"; 373 375 } 374 376 } 375 377 376 378 # basic ensemble stats 377 my $mean 378 my $meanStdev 379 my $stdev 380 my $binStdev 379 my $mean = $stats->value_for_flag ("-bg"); 380 my $meanStdev = $stats->value_for_flag ("-bg_mean_stdev"); 381 my $stdev = $stats->value_for_flag ("-bg_stdev"); 382 my $binStdev = $stats->value_for_flag ("-bin_stdev"); 381 383 my $fringe_mean = $stats->value_for_flag ("-fringe_0"); 382 384 my $fringe_err = $stats->value_for_flag ("-fringe_1"); … … 387 389 388 390 # other stats (flux depends on bg and exp_time) 389 my $fluxStats = Statistics::Descriptive::Sparse->new(); 391 my $fluxStats = Statistics::Descriptive::Sparse->new(); # Statistics calculator for means 390 392 $fluxStats->add_data(@fluxes); 391 my $flux = $fluxStats->mean(); 393 my $flux = $fluxStats->mean(); # Mean of the imfile means 392 394 393 395 # other stats … … 406 408 # for a bias). in general, use with ADDITIVE components 407 409 if ($reject_exp_mean > 0) { 408 409 410 411 412 413 414 415 } else { 416 410 if (abs($mean) > $reject_exp_mean * $stdev) { 411 print $logFile "Rejecting exposure based on bad mean counts: "; 412 $reject = 1; 413 } else { 414 print $logFile "Exposure mean counts meets requirements: "; 415 } 416 print $logFile "mean: $mean, stdev: $stdev (limit is: " . $reject_exp_mean * $stdev . ") \n"; 417 } else { 418 print $logFile "No rejection for exp mean\n"; 417 419 } 418 420 … … 423 425 # for a bias). in general, use with ADDITIVE components 424 426 if ($reject_exp_flux > 0) { 425 426 427 428 429 430 431 432 } else { 433 427 if (abs($flux) > $reject_exp_flux * $stdev) { 428 print $logFile "Rejecting exposure based on bad mean flux: "; 429 $reject = 1; 430 } else { 431 print $logFile "Exposure mean flux meets requirements: "; 432 } 433 print $logFile "flux: $flux, stdev: $stdev (limit is: " . $reject_exp_flux * $stdev . ")\n"; 434 } else { 435 print $logFile "No rejection for exp mean\n"; 434 436 } 435 437 … … 437 439 # limit 438 440 if ($reject_exp_stdev > 0) { 439 440 441 442 443 444 445 446 } else { 447 441 if ($stdev > $reject_exp_stdev) { 442 print $logFile "Rejecting exposure based on bad stdev: "; 443 $reject = 1; 444 } else { 445 print $logFile "Exposure stdev meets requirements: "; 446 } 447 print $logFile "$stdev vs $reject_exp_stdev\n"; 448 } else { 449 print $logFile "No rejection for exp stdev\n"; 448 450 } 449 451 … … 451 453 # than the limit 452 454 if ($reject_exp_meanstdev > 0) { 453 454 455 456 457 458 459 460 } else { 461 455 if ($meanStdev > $reject_exp_meanstdev) { 456 print $logFile "Rejecting exposure based on bad mean stdev: "; 457 $reject = 1; 458 } else { 459 print $logFile "Exposure mean stdev meets requirements: "; 460 } 461 print $logFile "$meanStdev vs $reject_exp_meanstdev\n"; 462 } else { 463 print $logFile "No rejection for exp mean stdev\n"; 462 464 } 463 465 464 466 # reject if the signal-to-noise is insufficient 465 467 if ($reject_exp_snr > 0) { 466 467 468 469 470 471 472 473 } else { 474 468 if (abs($mean) < abs($stdev * $reject_exp_snr)) { 469 print $logFile "Rejecting exposure based on poor S/N: \n"; 470 $reject = 1; 471 } else { 472 print $logFile "Exposure S/N meets requirements: \n"; 473 } 474 print $logFile "signal: $mean vs noise: $stdev (s/n limit is: $reject_exp_snr)\n"; 475 } else { 476 print $logFile "No rejection for exp S/N\n"; 475 477 } 476 478 # reject if the exposure ensemble stdev is deviant 477 479 if ($reject_exp_bin_stdev > 0) { 478 479 480 481 482 483 484 485 } else { 486 480 if ($binStdev > $reject_exp_bin_stdev) { 481 print $logFile "Rejecting exposure based on bad binned stdev: "; 482 $reject = 1; 483 } else { 484 print $logFile "Exposure binned stdev meets requirements: "; 485 } 486 print $logFile "$stdev vs $reject_exp_bin_stdev\n"; 487 } else { 488 print $logFile "No rejection for exp binned stdev\n"; 487 489 } 488 490 # reject if the signal-to-noise is insufficient 489 491 if ($reject_exp_bin_snr > 0) { 490 491 492 493 494 495 496 497 } else { 498 492 if (abs($mean) < abs($binStdev * $reject_exp_bin_snr)) { 493 print $logFile "Rejecting exposure based on poor binned S/N: \n"; 494 $reject = 1; 495 } else { 496 print $logFile "Exposure binned S/N meets requirements: \n"; 497 } 498 print $logFile "signal: $mean vs noise: $binStdev (s/n limit is: $reject_exp_bin_snr)\n"; 499 } else { 500 print $logFile "No rejection for exp binned S/N\n"; 499 501 } 500 502 … … 510 512 unless ($no_update) { 511 513 ( $success, $error_code, $full_buf, $stdout_buf, $stderr_buf ) = 512 514 run(command => $command, verbose => $verbose); 513 515 unless ($success) { 514 515 516 516 $error_code = (($error_code >> 8) or $PS_EXIT_PROG_ERROR); 517 warn("Unable to perform dettool -addresidexp: $error_code\n"); 518 exit($error_code); 517 519 } 518 520 } else { … … 531 533 { 532 534 my $msg = shift; # Warning message on die 533 my $det_id = shift; 534 my $iter = shift; 535 my $det_id = shift; # Detrend identifier 536 my $iter = shift; # Iteration 535 537 my $exp_id = shift; # Exposure tag 536 538 my $exit_code = shift; # Exit code to add … … 538 540 carp($msg); 539 541 if (defined $det_id and defined $iter and defined $exp_id and not $no_update) { 540 541 542 543 544 545 542 my $command = "$dettool -addresidexp"; 543 $command .= " -det_id $det_id"; 544 $command .= " -iteration $iter"; 545 $command .= " -exp_id $exp_id"; 546 $command .= " -code $exit_code"; 547 $command .= " -dbname $dbname" if defined $dbname; 546 548 system ($command); 547 549 } … … 552 554 sub rejection_limit 553 555 { 554 my $name = shift; # Rejection limit to555 my $type = shift; 556 my $filter = shift; 556 my $name = shift; # Rejection limit to 557 my $type = shift; # Type of exposure 558 my $filter = shift; # Filter 557 559 558 560 my $value = $ipprc->rejection( $name, $det_type, $filter ); 559 561 if (not defined $value) { 560 561 562 $filter = "(no filter)" if not defined $filter; 563 die "Unable to determine $name rejection limit for $det_type with $filter.\n"; 562 564 } 563 565 … … 589 591 ## 5) is the mean stdev > limit (poor images) 590 592 ## 6) calculate the ensemble stdev (rms) of remaining (excluding rejects) 591 ## 7) is the ensemble stdev > limit 593 ## 7) is the ensemble stdev > limit 592 594 593 595 ## for flats, the value used for stdev should be the fractional stdev
Note:
See TracChangeset
for help on using the changeset viewer.
