Changeset 36441

Timestamp:

Jan 14, 2014, 2:28:47 PM (12 years ago)

Author:

bills

Message:

Iniital implementation of the full force stages to the pipeline.
Maginitude and galactic coordinate limits on extended source fits

Location:

trunk

Files:

• . (modified) (1 prop)
• Nebulous (modified) (1 prop)
• Nebulous-Server (modified) (1 prop)
• Ohana (modified) (1 prop)
• Ohana/src/dvomerge (modified) (1 prop)
• Ohana/src/libohana/src (modified) (1 prop)
• Ohana/src/opihi (modified) (1 prop)
• Ohana/src/opihi/cmd.astro (modified) (1 prop)
• Ohana/src/opihi/cmd.data (modified) (1 prop)
• Ohana/src/relastro/src (modified) (1 prop)
• dbconfig/cam.md (modified) (1 prop)
• dbconfig/changes.txt (modified) (1 diff, 1 prop)
• dbconfig/ff.md (copied) (copied from branches/bills_branches/bills_201312/dbconfig/ff.md )
• dbconfig/ipp.m4 (modified) (1 diff)
• ippMonitor (modified) (1 prop)
• ippScripts/Build.PL (modified) (1 diff)
• ippScripts/scripts/chip_imfile.pl (modified) (1 prop)
• ippScripts/scripts/destreak_restore_camera.pl (modified) (1 prop)
• ippScripts/scripts/ipp_apply_burntool_single.pl (modified) (1 prop)
• ippScripts/scripts/lap_science.pl (modified) (1 prop)
• ippScripts/scripts/magic_destreak.pl (modified) (1 prop)
• ippScripts/scripts/psphot_fullforce_summary.pl (copied) (copied from branches/bills_branches/bills_201312/ippScripts/scripts/psphot_fullforce_summary.pl )
• ippScripts/scripts/psphot_fullforce_warp.pl (modified) (18 diffs)
• ippScripts/scripts/publish_file.pl (modified) (1 prop)
• ippScripts/scripts/staticsky.pl (modified) (4 diffs)
• ippScripts/scripts/warp_skycell.pl (modified) (1 prop)
• ippTasks (modified) (1 prop)
• ippTasks/Makefile.am (modified) (1 diff)
• ippTasks/fullforce.pro (copied) (copied from branches/bills_branches/bills_201312/ippTasks/fullforce.pro )
• ippToPsps (modified) (1 prop)
• ippTools/share (modified) (1 prop)
• ippTools/share/Makefile.am (modified) (1 diff)
• ippTools/share/camtool_find_pendingimfile.sql (modified) (1 prop)
• ippTools/share/chiptool_setimfiletoupdate.sql (modified) (1 prop)
• ippTools/share/fftool_definebyquery.sql (copied) (copied from branches/bills_branches/bills_201312/ippTools/share/fftool_definebyquery.sql )
• ippTools/share/fftool_definebyquery_select_warps.sql (copied) (copied from branches/bills_branches/bills_201312/ippTools/share/fftool_definebyquery_select_warps.sql )
• ippTools/share/fftool_result.sql (copied) (copied from branches/bills_branches/bills_201312/ippTools/share/fftool_result.sql )
• ippTools/share/fftool_revert.sql (copied) (copied from branches/bills_branches/bills_201312/ippTools/share/fftool_revert.sql )
• ippTools/share/fftool_revertsummary.sql (copied) (copied from branches/bills_branches/bills_201312/ippTools/share/fftool_revertsummary.sql )
• ippTools/share/fftool_summary.sql (copied) (copied from branches/bills_branches/bills_201312/ippTools/share/fftool_summary.sql )
• ippTools/share/fftool_summaryinputs.sql (copied) (copied from branches/bills_branches/bills_201312/ippTools/share/fftool_summaryinputs.sql )
• ippTools/share/fftool_toadvance.sql (copied) (copied from branches/bills_branches/bills_201312/ippTools/share/fftool_toadvance.sql )
• ippTools/share/fftool_todo.sql (copied) (copied from branches/bills_branches/bills_201312/ippTools/share/fftool_todo.sql )
• ippTools/share/pxadmin_create_tables.sql (modified) (1 diff, 1 prop)
• ippTools/share/pxadmin_drop_tables.sql (modified) (1 diff)
• ippTools/share/warptool_towarped.sql (modified) (1 prop)
• ippTools/src (modified) (2 props)
• ippTools/src/Makefile.am (modified) (3 diffs)
• ippTools/src/camtool.c (modified) (1 prop)
• ippTools/src/camtoolConfig.c (modified) (1 prop)
• ippTools/src/fftool.c (copied) (copied from branches/bills_branches/bills_201312/ippTools/src/fftool.c )
• ippTools/src/fftool.h (copied) (copied from branches/bills_branches/bills_201312/ippTools/src/fftool.h )
• ippTools/src/fftoolConfig.c (copied) (copied from branches/bills_branches/bills_201312/ippTools/src/fftoolConfig.c )
• ippTools/src/magictool.c (modified) (1 prop)
• ippconfig (modified) (1 prop)
• ippconfig/gpc1 (modified) (1 prop)
• ippconfig/recipes/filerules-split.mdc (modified) (1 diff)
• ippconfig/recipes/ppSub.config (modified) (1 prop)
• ippconfig/recipes/psastro.config (modified) (2 diffs)
• ippconfig/recipes/psphot.config (modified) (6 diffs)
• ippconfig/recipes/reductionClasses.mdc (modified) (1 prop)
• ppImage/src (modified) (1 prop)
• ppSub/src/ppSubMatchPSFs.c (modified) (1 diff)
• psLib/src/math/psMixtureModels.c (modified) (1 prop)
• psModules (modified) (1 prop)
• psModules/src/camera/pmFPAfile.c (modified) (1 diff)
• psModules/src/objects/pmModelClass.c (modified) (3 diffs)
• psModules/src/objects/pmModelClass.h (modified) (1 diff)
• psModules/src/objects/pmSourceIO.c (modified) (7 diffs)
• psModules/src/objects/pmSourceIO.h (modified) (1 diff)
• psModules/src/objects/pmSourceIO_CFF.c (modified) (8 diffs)
• psModules/src/objects/pmSourceIO_CMF.c.in (modified) (9 diffs)
• psModules/src/objects/pmSourceIO_PS1_CAL_0.c (modified) (1 diff)
• psModules/src/objects/pmSourceIO_PS1_DEV_0.c (modified) (1 diff)
• psModules/src/objects/pmSourceIO_PS1_DEV_1.c (modified) (1 diff)
• psModules/src/objects/pmSourceIO_SMPDATA.c (modified) (1 diff)
• psastro/src/psastroDataSave.c (modified) (2 diffs)
• psastro/src/psastroDefineFiles.c (modified) (1 diff)
• psconfig (modified) (1 prop)
• psphot (modified) (1 prop)
• psphot/src (modified) (2 props)
• psphot/src/Makefile.am (modified) (3 diffs)
• psphot/src/psmakecff.c (modified) (8 diffs)
• psphot/src/psphot.h (modified) (1 diff)
• psphot/src/psphotChooseAnalysisOptions.c (modified) (19 diffs)
• psphot/src/psphotFullForceSummary.c (copied) (copied from branches/bills_branches/bills_201312/psphot/src/psphotFullForceSummary.c )
• psphot/src/psphotFullForceSummaryReadout.c (copied) (copied from branches/bills_branches/bills_201312/psphot/src/psphotFullForceSummaryReadout.c )
• psphot/src/psphotGalaxyShape.c (modified) (3 diffs)
• psphot/src/psphotStackImageLoop.c (modified) (1 prop)
• pstamp/scripts (modified) (1 prop)
• psvideophot (modified) (1 prop)

trunk/dbconfig/changes.txt

@@ -2432 +2432 @@
 UPDATE dbversion set schema_version = '1.1.77', updated= CURRENT_TIMESTAMP();
 
+--
+
+CREATE TABLE fullForceRun (
+    ff_id           BIGINT NOT NULL AUTO_INCREMENT,
+    skycal_id       BIGINT,
+    sources_path_base VARCHAR(255),
+    state           VARCHAR(64),
+    workdir         VARCHAR(255),
+    label           VARCHAR(64),
+    data_group      VARCHAR(64),
+    dist_group      VARCHAR(64),
+    note            VARCHAR(255),
+    reduction       VARCHAR(64),
+    registered TIMESTAMP DEFAULT CURRENT_TIMESTAMP, -- time run was registered
+    PRIMARY KEY(ff_id),
+    KEY(state),
+    KEY(label),
+    KEY(data_group),
+    KEY(skycal_id),
+    FOREIGN KEY(skycal_id) REFERENCES skycalRun(skycal_id)
+) ENGINE=innodb DEFAULT CHARSET=latin1;
+
+CREATE TABLE fullForceInput (
+    ff_id           BIGINT,
+    warp_id         BIGINT,
+    PRIMARY KEY(ff_id, warp_id),
+    FOREIGN KEY(ff_id) REFERENCES fullForceRun(ff_id),
+    FOREIGN KEY(warp_id) REFERENCES warpRun(warp_id)
+) ENGINE=innodb DEFAULT CHARSET=latin1;
+
+CREATE TABLE fullForceResult (
+    ff_id           BIGINT,
+    warp_id         BIGINT,
+    path_base       VARCHAR(255) NOT NULL,
+    dtime_script    FLOAT,
+    quality         SMALLINT NOT NULL,
+    hostname        VARCHAR(64) NOT NULL,
+    software_ver    VARCHAR(16),
+    fault           SMALLINT NOT NULL,
+    PRIMARY KEY(ff_id, warp_id),
+    KEY(fault),
+    KEY(quality),
+    FOREIGN KEY(ff_id) REFERENCES fullForceRun(ff_id),
+    FOREIGN KEY(warp_id) REFERENCES warpRun(warp_id)
+) ENGINE=innodb DEFAULT CHARSET=latin1;
+
+CREATE TABLE fullForceSummary (
+    ff_id           BIGINT,
+    path_base       VARCHAR(255) NOT NULL,
+    dtime_script    FLOAT,
+    quality         SMALLINT NOT NULL,
+    hostname        VARCHAR(64) NOT NULL,
+    software_ver    VARCHAR(16),
+    fault           SMALLINT NOT NULL,
+    PRIMARY KEY(ff_id),
+    KEY(fault),
+    KEY(quality),
+    FOREIGN KEY(ff_id) REFERENCES fullForceRun(ff_id)
+) ENGINE=innodb DEFAULT CHARSET=latin1;
+
+UPDATE dbversion set schema_version = '1.1.78', updated= CURRENT_TIMESTAMP();
+

trunk/dbconfig/ipp.m4

@@ -40 +40 @@
 include(skycell.md)
 include(release.md)
+include(ff.md)

trunk/ippScripts/Build.PL

@@ -133 +133 @@
         scripts/queuestaticsky.pl
         scripts/psphot_fullforce_warp.pl
+        scripts/psphot_fullforce_summary.pl
     )],
     dist_abstract => 'Scripts for running the Pan-STARRS IPP',

trunk/ippScripts/scripts/psphot_fullforce_warp.pl

@@ -31 +31 @@
 my $psphotFullForce = can_run('psphotFullForce') or (warn "Can't find psphotFullForce" and $missing_tools = 1);
 my $ppStatsFromMetadata = can_run('ppStatsFromMetadata') or (warn "Can't find ppStatsFromMetadata" and $missing_tools = 1);
-# XXX: fftool is yet to be written
-my $fftool = "fftool";  # can_run('fftool') or (warn "Can't find fftool" and $missing_tools = 1);
+my $fftool = can_run('fftool') or (warn "Can't find fftool" and $missing_tools = 1);
 if ($missing_tools) {
     warn("Can't find required tools.");
@@ -38 +37 @@
 }
 
-my ($ffw_id, $warp_id, $skycell_id, $path_base, $sourceroot, $camera);
+my ($ff_id, $warp_id, $skycell_id, $path_base, $sourceroot, $camera);
 my ($outroot, $reduction);
 my ($dbname, $threads, $verbose, $no_update, $no_op, $redirect);
 
 GetOptions(
-    'ffw_id=s'          => \$ffw_id,
+    'ff_id=s'          => \$ff_id,
     'warp_id=s'         => \$warp_id,   # warp identifier
     'skycell_id=s'      => \$skycell_id,# Skycell identifier
     'warp_path_base=s'  => \$path_base, # path_base of the warp skycell
-    'sourceroot=s'      => \$sourceroot,# path_base of sources
+    'sources_path_base=s' => \$sourceroot,# path_base of sources
     'camera=s'          => \$camera,    # camera name of sources
     'dbname|d=s'        => \$dbname,    # Database name
@@ -61 +60 @@
 pod2usage( -msg => "Unknown option: @ARGV", -exitval => 2 ) if @ARGV;
 pod2usage(
-    -msg => "Required options: --ffw_id --sourceroot --skycell_id --warp_path_base --outroot --camera",
+    -msg => "Required options: --ff_id --warp_id --sourceroot --skycell_id --warp_path_base --outroot --camera",
     -exitval => 3,
-          ) unless defined $ffw_id,
+          ) unless defined $ff_id,
     and defined $sourceroot
-    and (defined $path_base or defined $warp_id) # if we don't have warp's path_base we need warp_id
+    and defined $path_base 
+    and defined $warp_id
     and defined $skycell_id
     and defined $camera
     and defined $outroot;
 
-# XXX: fftool is not ready to run commands that update the database
-$no_update = 1;
-
-my $ipprc = PS::IPP::Config->new($camera) or my_die( "Unable to set up", $ffw_id, $skycell_id, $PS_EXIT_CONFIG_ERROR );
+my $ipprc = PS::IPP::Config->new($camera) or my_die( "Unable to set up", $ff_id, $warp_id, $skycell_id, $PS_EXIT_CONFIG_ERROR );
 
 my $neb;
@@ -84 +81 @@
 
 $ipprc->redirect_to_logfile($logDest) or my_die( "Unable to redirect output", 
-    $ffw_id, $skycell_id, $PS_EXIT_SYS_ERROR ) if $redirect;
-
-if (!$path_base) {
+    $ff_id, $warp_id, $skycell_id, $PS_EXIT_SYS_ERROR ) if $redirect;
+
+if (0) { 
+# if (!$path_base) {
     # If path_base is not supplied, look it up in the database.
     my $mdcParser = PS::IPP::Metadata::Config->new; # Parser for metadata config files
@@ -98 +96 @@
             $error_code = (($error_code >> 8) or $PS_EXIT_PROG_ERROR);
             &my_die("Unable to perform warptool -warpskyfile -inputskyfile: $error_code", 
-                $ffw_id, $skycell_id, $error_code);
+                $ff_id, $warp_id, $skycell_id, $error_code);
         }
 
         my $metadata = $mdcParser->parse(join "", @$stdout_buf) or
-            &my_die("Unable to parse metadata config doc", $ffw_id, $skycell_id, $PS_EXIT_PROG_ERROR);
+            &my_die("Unable to parse metadata config doc", $ff_id, $warp_id, $skycell_id, $PS_EXIT_PROG_ERROR);
         $files = parse_md_list($metadata) or
-            &my_die("Unable to parse metadata list", $ffw_id, $skycell_id, $PS_EXIT_PROG_ERROR);
-    }
-
-    &my_die("Input list does not contain exactly one elements", $ffw_id, $skycell_id, $PS_EXIT_SYS_ERROR) 
+            &my_die("Unable to parse metadata list", $ff_id, $warp_id, $skycell_id, $PS_EXIT_PROG_ERROR);
+    }
+
+    &my_die("Input list does not contain exactly one elements", $ff_id, $warp_id, $skycell_id, $PS_EXIT_SYS_ERROR) 
         unless scalar @$files == 1;
 
@@ -113 +111 @@
 
     $path_base = $warp->{path_base};
-    &my_die("Couldn't find input path in warptool output", $ffw_id, $skycell_id, $PS_EXIT_UNKNOWN_ERROR) 
+    &my_die("Couldn't find input path in warptool output", $ff_id, $warp_id, $skycell_id, $PS_EXIT_UNKNOWN_ERROR) 
         unless defined $path_base;
 
@@ -123 +121 @@
 unless ($recipe_psphot) {
     &my_die("Couldn't find selected reduction for PSPHOT: $reduction\n", 
-        $ffw_id, $skycell_id, $PS_EXIT_CONFIG_ERROR);
+        $ff_id, $warp_id, $skycell_id, $PS_EXIT_CONFIG_ERROR);
 }
 
@@ -134 +132 @@
 print "recipe_psphot: $recipe_psphot\n";
 
+# use psf measured on input warp
+# XXX: get this from recipe
+my $useWarpPSF = 0;
+
 my $input         = $ipprc->filename('PSWARP.OUTPUT', $path_base);
 my $inputMask     = $ipprc->filename('PSWARP.OUTPUT.MASK', $path_base);
 my $inputVariance = $ipprc->filename('PSWARP.OUTPUT.VARIANCE', $path_base);
-my $inputPSF      = $ipprc->filename('PSPHOT.PSF.SKY.SAVE', $path_base);
+my $inputPSF      = $useWarpPSF ? $ipprc->filename('PSPHOT.PSF.SKY.SAVE', $path_base) : "";
 my $inputSources  = $ipprc->filename('PSPHOT.OUTPUT.CFF', $sourceroot);
 
@@ -144 +146 @@
     print "inputMask:     $inputMask\n";
     print "inputVariance: $inputVariance\n";
-    # print "inputPath:     $path_base\n";
-    print "inputPSF:      $inputPSF\n";
+    print "inputPSF:      $inputPSF\n" if $inputPSF;
     print "inputSources:  $inputSources\n";
 }
 
 # check that the inputs exist (and have non-zero size)
-&my_die("Couldn't find input: $input", $ffw_id, $skycell_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists($input);
-&my_die("Couldn't find input: $inputMask", $ffw_id, $skycell_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists($inputMask);
-&my_die("Couldn't find input: $inputVariance", $ffw_id, $skycell_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists($inputVariance);
-&my_die("Couldn't find input: $inputPSF", $ffw_id, $skycell_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists($inputPSF);
-&my_die("Couldn't find input: $inputSources", $ffw_id, $skycell_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists($inputSources);
+&my_die("Couldn't find input: $input", $ff_id, $warp_id, $skycell_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists($input);
+&my_die("Couldn't find input: $inputMask", $ff_id, $warp_id, $skycell_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists($inputMask);
+&my_die("Couldn't find input: $inputVariance", $ff_id, $warp_id, $skycell_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists($inputVariance);
+&my_die("Couldn't find input: $inputPSF", $ff_id, $warp_id, $skycell_id, $PS_EXIT_SYS_ERROR) if ($inputPSF && !$ipprc->file_exists($inputPSF));
+&my_die("Couldn't find input: $inputSources", $ff_id, $warp_id, $skycell_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists($inputSources);
 
 my $dump_config = 1; 
@@ -193 +194 @@
         unless ($success) {
             $error_code = (($error_code >> 8) or $PS_EXIT_PROG_ERROR);
-            &my_die("Unable to perform ppSub: $error_code", $ffw_id, $skycell_id, $error_code);
+            &my_die("Unable to perform ppSub: $error_code", $ff_id, $warp_id, $skycell_id, $error_code);
         }
 
@@ -207 +208 @@
             unless ($success) {
                 $error_code = (($error_code >> 8) or $PS_EXIT_PROG_ERROR);
-                &my_die("Unable to perform ppStatsFromMetadata: $error_code", $ffw_id, $skycell_id, $error_code);
+                &my_die("Unable to perform ppStatsFromMetadata: $error_code", $ff_id, $warp_id, $skycell_id, $error_code);
             }
             foreach my $line (@$stdout_buf) {
@@ -226 +227 @@
 # Add the result to the database
 {
-    my $command = "$fftool -ffw_id $ffw_id -skycell_id $skycell_id";
-    $command .= " -addwarped -path_base $outroot";
+    my $command = "$fftool -ff_id $ff_id -warp_id $warp_id"; 
+    $command .= " -addresult -path_base $outroot";
     $command .= " $cmdflags";
     $command .= (" -dtime_script " . ((DateTime->now->mjd - $mjd_start) * 86400));
@@ -239 +240 @@
             $error_code = (($error_code >> 8) or $PS_EXIT_PROG_ERROR);
             my $err_message = "Unable to perform fftool -addwarped" ;
-                &my_die("$err_message: $error_code", $ffw_id, $skycell_id, $error_code);
+                &my_die("$err_message: $error_code", $ff_id, $warp_id, $skycell_id, $error_code);
         }
     } else {
@@ -262 +263 @@
     my $error;
     my $output = $ipprc->prepare_output($filerule, $outroot, undef, $delete, \$error)
-                    or &my_die("failed to prepare output file for: $filerule", $ffw_id, $skycell_id, $error);
+                    or &my_die("failed to prepare output file for: $filerule", $ff_id, $warp_id, $skycell_id, $error);
     return $output;
 }
@@ -275 +276 @@
     }
 
-    &my_die("Couldn't find expected output file: $file",  $ffw_id, $skycell_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists($file);
+    &my_die("Couldn't find expected output file: $file",  $ff_id, $warp_id, $skycell_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists($file);
 
     # Funpack to confirm we've really made things correctly
     my $diskfile = $ipprc->file_resolve($file);
     if ($diskfile =~ /fits/) {
-        my $funpack  = can_run('funpack') or &my_die ("Can't find funpack", $ffw_id, $skycell_id, $PS_EXIT_SYS_ERROR);
+        my $funpack  = can_run('funpack') or &my_die ("Can't find funpack", $ff_id, $warp_id, $skycell_id, $PS_EXIT_SYS_ERROR);
         my $check_command = "$funpack -S $diskfile > /dev/null";
         my ( $success, $error_code, $full_buf, $stdout_buf, $stderr_buf ) =
             run(command => $check_command, verbose => $verbose);
         if (!$success) {
-            &my_die("Output file not a valid fits file: $file", $ffw_id, $skycell_id, $PS_EXIT_SYS_ERROR);
+            &my_die("Output file not a valid fits file: $file", $ff_id, $warp_id, $skycell_id, $PS_EXIT_SYS_ERROR);
         }
     }
@@ -291 +292 @@
 
     if ($replicate and $neb) {
-        $ipprc->replicate_file($file) or &my_die("failed to replicate: $file\n",  $ffw_id, $skycell_id, $PS_EXIT_SYS_ERROR);
+        $ipprc->replicate_file($file) or &my_die("failed to replicate: $file\n",  $ff_id, $warp_id, $skycell_id, $PS_EXIT_SYS_ERROR);
     }
 }
@@ -299 +300 @@
 {
     my $msg = shift;            # Warning message on die
-    my $ffw_id = shift;         # full force warp identifier
+    my $ff_id = shift;          # full force run identifier
+    my $warp_id = shift;        # full force warp id
     my $skycell_id = shift;     # Skycell identifier
     my $exit_code = shift;      # Exit code to add
@@ -306 +308 @@
 
     warn($msg);
-    if (defined $ffw_id and defined $skycell_id) {
-        my $command = "$fftool -ffw_id $ffw_id -skycell_id $skycell_id -fault $exit_code";
-        $command .= " -addffskyfile";
+    if (defined $ff_id and defined $skycell_id) {
+        my $command = "$fftool -ff_id $ff_id -warp_id $warp_id -fault $exit_code";
+        $command .= " -addresult";
         $command .= (" -dtime_script " . ((DateTime->now->mjd - $mjd_start) * 86400));
         $command .= " -hostname $host" if defined $host;

trunk/ippScripts/scripts/staticsky.pl

@@ -138 +138 @@
         my $configuration = $ipprc->filename("PSPHOT.STACK.CONFIG", $outroot);
 
+        my $needConvolvedImages = 0;
+
         foreach my $file (@$files) {
             print $listFile "INPUT   METADATA\n";
@@ -146 +148 @@
             my $stack_id = $file->{stack_id};
 
-            my $imageCnv  = $ipprc->filename("PPSTACK.OUTPUT",          $path_base ); # Image name
-            my $maskCnv   = $ipprc->filename("PPSTACK.OUTPUT.MASK",     $path_base ); # Mask name
-            my $weightCnv = $ipprc->filename("PPSTACK.OUTPUT.VARIANCE", $path_base ); # Weight name
-            my $expnumCnv = $ipprc->filename("PPSTACK.OUTPUT.EXPNUM",   $path_base ); # Expnum name
+            my ($imageCnv, $maskCnv, $weightCnv, $expnumCnv);
+            if ($needConvolvedImages) {
+                $imageCnv  = $ipprc->filename("PPSTACK.OUTPUT",          $path_base ); # Image name
+                $maskCnv   = $ipprc->filename("PPSTACK.OUTPUT.MASK",     $path_base ); # Mask name
+                $weightCnv = $ipprc->filename("PPSTACK.OUTPUT.VARIANCE", $path_base ); # Weight name
+                $expnumCnv = $ipprc->filename("PPSTACK.OUTPUT.EXPNUM",   $path_base ); # Expnum name
+            }
 
             my $imageRaw  = $ipprc->filename("PPSTACK.UNCONV",          $path_base ); # Image name
@@ -166 +171 @@
             &my_die("Couldn't find input: $weightRaw", $sky_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists("$weightRaw");
             &my_die("Couldn't find input: $expnumRaw", $sky_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists("$expnumRaw");
-            &my_die("Couldn't find input: $imageCnv",  $sky_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists("$imageCnv");
-            &my_die("Couldn't find input: $maskCnv",   $sky_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists("$maskCnv");
-            &my_die("Couldn't find input: $weightCnv", $sky_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists("$weightCnv");
-            &my_die("Couldn't find input: $expnumCnv", $sky_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists("$expnumCnv");
-            &my_die("Couldn't find input: $psfCnv",    $sky_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists("$psfCnv");
+            if ($needConvolvedImages) {
+                &my_die("Couldn't find input: $imageCnv",  $sky_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists("$imageCnv");
+                &my_die("Couldn't find input: $maskCnv",   $sky_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists("$maskCnv");
+                &my_die("Couldn't find input: $weightCnv", $sky_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists("$weightCnv");
+                &my_die("Couldn't find input: $expnumCnv", $sky_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists("$expnumCnv");
+                &my_die("Couldn't find input: $psfCnv",    $sky_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists("$psfCnv");
+            }
             &my_die("Couldn't find input: $sources",   $sky_id, $PS_EXIT_SYS_ERROR) unless $ipprc->file_exists("$sources");
 
@@ -178 +185 @@
             print $listFile "  RAW:VARIANCE  STR  " . $weightRaw . "\n";
             print $listFile "  RAW:EXPNUM    STR  " . $expnumRaw . "\n";
-
-            print $listFile "  CNV:IMAGE     STR  " . $imageCnv  . "\n";
-            print $listFile "  CNV:MASK      STR  " . $maskCnv   . "\n";
-            print $listFile "  CNV:VARIANCE  STR  " . $weightCnv . "\n";
-            print $listFile "  CNV:EXPNUM    STR  " . $expnumCnv . "\n";
-            print $listFile "  CNV:PSF       STR  " . $psfCnv    . "\n";
-
             print $listFile "  SOURCES       STR  " . $sources   . "\n";
+            if ($needConvolvedImages) {
+                print $listFile "  CNV:IMAGE     STR  " . $imageCnv  . "\n";
+                print $listFile "  CNV:MASK      STR  " . $maskCnv   . "\n";
+                print $listFile "  CNV:VARIANCE  STR  " . $weightCnv . "\n";
+                print $listFile "  CNV:EXPNUM    STR  " . $expnumCnv . "\n";
+                print $listFile "  CNV:PSF       STR  " . $psfCnv    . "\n";
+            }
 
             print $listFile "END\n\n";

trunk/ippTasks/Makefile.am

@@ -47 +47 @@
         lap.pro \
         vp.pro \
-        bg.regeneration.pro
+        bg.regeneration.pro \
+        fullforce.pro
 
 other_files = \

trunk/ippTools/share/Makefile.am

@@ -488 +488 @@
         releasetool_definerelgroup_select_exp_lap.sql \
         releasetool_stacksummary.sql \
-        releasetool_pendingrelgroup.sql
+        releasetool_pendingrelgroup.sql \
+        fftool_definebyquery.sql \
+        fftool_definebyquery_select_warps.sql \
+        fftool_result.sql \
+        fftool_revert.sql \
+        fftool_revertsummary.sql \
+        fftool_summary.sql \
+        fftool_todo.sql \
+        fftool_toadvance.sql
 

trunk/ippTools/share/pxadmin_create_tables.sql

@@ -2307 +2307 @@
 ) ENGINE=InnoDB DEFAULT CHARSET=latin1;
 
+CREATE TABLE fullForceRun (
+    ff_id           BIGINT NOT NULL AUTO_INCREMENT,
+    skycal_id       BIGINT,
+    sources_path_base VARCHAR(255),
+    state           VARCHAR(64),
+    workdir         VARCHAR(255),
+    label           VARCHAR(64),
+    data_group      VARCHAR(64),
+    dist_group      VARCHAR(64),
+    note            VARCHAR(255),
+    reduction       VARCHAR(64),
+    registered TIMESTAMP DEFAULT CURRENT_TIMESTAMP, -- time run was registered
+    PRIMARY KEY(ff_id),
+    KEY(state),
+    KEY(label),
+    KEY(data_group),
+    KEY(skycal_id),
+    FOREIGN KEY(skycal_id) REFERENCES skycalRun(skycal_id)
+) ENGINE=innodb DEFAULT CHARSET=latin1;
+
+CREATE TABLE fullForceInput (
+    ff_id           BIGINT,
+    warp_id         BIGINT,
+    PRIMARY KEY(ff_id, warp_id),
+    FOREIGN KEY(ff_id) REFERENCES fullForceRun(ff_id),
+    FOREIGN KEY(warp_id) REFERENCES warpRun(warp_id)
+) ENGINE=innodb DEFAULT CHARSET=latin1;
+
+CREATE TABLE fullForceResult (
+    ff_id           BIGINT,
+    warp_id         BIGINT,
+    path_base       VARCHAR(255) NOT NULL,
+    dtime_script    FLOAT,
+    quality         SMALLINT NOT NULL,
+    hostname        VARCHAR(64) NOT NULL,
+    software_ver    VARCHAR(16),
+    fault           SMALLINT NOT NULL,
+    PRIMARY KEY(ff_id, warp_id),
+    KEY(fault),
+    KEY(quality),
+    FOREIGN KEY(ff_id) REFERENCES fullForceRun(ff_id),
+    FOREIGN KEY(warp_id) REFERENCES warpRun(warp_id)
+) ENGINE=innodb DEFAULT CHARSET=latin1;
+
+CREATE TABLE fullForceSummary (
+    ff_id           BIGINT,
+    path_base       VARCHAR(255) NOT NULL,
+    dtime_script    FLOAT,
+    quality         SMALLINT NOT NULL,
+    hostname        VARCHAR(64) NOT NULL,
+    software_ver    VARCHAR(16),
+    fault           SMALLINT NOT NULL,
+    PRIMARY KEY(ff_id),
+    KEY(fault),
+    KEY(quality),
+    FOREIGN KEY(ff_id) REFERENCES fullForceRun(ff_id)
+) ENGINE=innodb DEFAULT CHARSET=latin1;
+
 
 -- These comment lines are here to avoid an empty query error.

trunk/ippTools/share/pxadmin_drop_tables.sql

@@ -120 +120 @@
 DROP TABLE IF EXISTS dqstatsContent;
 DROP TABLE IF EXISTS dqstatsRun;
+DROP TABLE IF EXISTS fullForceRun;
+DROP TABLE IF EXISTS fullForceInput;
+DROP TABLE IF EXISTS fullForceResult;
+DROP TABLE IF EXISTS fullForceSummary;
+
 
 SET FOREIGN_KEY_CHECKS=1

trunk/ippTools/src

@@ -47 +47 @@
 sctool
 releasetool
+fftool

@@ -47 +47 @@
 sctool
 releasetool
+fftool

trunk/ippTools/src/Makefile.am

@@ -33 +33 @@
         vptool \
         sctool \
-        releasetool
+        releasetool \
+        fftool
 
 pkginclude_HEADERS = \
@@ -87 +88 @@
         vptool.h \
         sctool.h \
-        releasetool.h
+        releasetool.h \
+        fftool.h
 
 lib_LTLIBRARIES = libpxtools.la
@@ -332 +334 @@
     releasetoolConfig.c
 
+fftool_CFLAGS = $(PSLIB_CFLAGS) $(PSMODULES_CFLAGS) $(IPPDB_CFLAGS)
+fftool_LDADD = $(PSLIB_LIBS) $(PSMODULES_LIBS) $(IPPDB_LIBS) libpxtools.la
+fftool_SOURCES = \
+    fftool.c \
+    fftoolConfig.c
+
 clean-local:
         -rm -f TAGS

trunk/ippconfig/recipes/filerules-split.mdc

@@ -262 +262 @@
 PSPHOT.STACK.RESID           OUTPUT {OUTPUT}.stk.{FILE.ID}.res.fits   IMAGE           COMP_SUB   FPA        TRUE      NONE
 PSPHOT.STACK.CONFIG          OUTPUT {OUTPUT}.psphotStack.mdc          TEXT            NONE       FPA        TRUE      NONE
+
+PSPHOT.FULLFORCE.OUTPUT      OUTPUT {OUTPUT}.cmf                      CMF             NONE       FPA        TRUE      NONE
                                                      
 SOURCE.PLOT.RAW.MOMENTS      OUTPUT {OUTPUT}.{CHIP.NAME}.mnt.png      KAPA            NONE       CHIP       TRUE      NONE

trunk/ippconfig/recipes/psastro.config

@@ -3 +3 @@
 PSASTRO.ONLY.REFSTARS      BOOL FALSE  # skip all but refstar matches
 PSASTRO.SAVE.REFMATCH      BOOL FALSE  # save refstar matches as table in output smf file
+PSASTRO.SAVE.CFF           BOOL FALSE  # save a cff file (input for psphotFullForce)
 
 # select which WCS style to use on output images.
@@ -279 +280 @@
 
 STATICSKY_CAL   METADATA
+    PSASTRO.SAVE.CFF           BOOL TRUE  # save a cff file (input for psphotFullForce)
 END
 

trunk/ippconfig/recipes/psphot.config

@@ -196 +196 @@
 EXT.NSIGMA.LIMIT.USE                BOOL  TRUE
 
-EXT.FIT.MIN.GAL.LIMIT               F32   10.0
-EXT.FIT.MIN.GAL.LIMIT.USE           BOOL  FALSE
+# Limits on extended source fits by galactic coordinates 
+EXT.FIT.MIN.GAL.LIMIT.USE           BOOL  FALSE # if true apply the galactic latitude cut
+# objects are skipped for extended source fits if they have a abs(galactic latitude) larger than
+# a limit value which is a function of galactic longitude.
+#     The galactic latitude limit is the sum of the following constant value ...
+EXT.FIT.MIN.GAL.LIMIT               F32   20.0
+# ... and a gaussian function of galactic longitude centered at b = 0 with this magnitude
+EXT.FIT.MIN.GAL.LIMIT.BULGE         F32   15.0  
+# ... and this sigma
+EXT.FIT.MIN.GAL.LIMIT.BULGE.SIGMA   F32   50.0  # with this sigma value
+
 
 KRON_ITERATIONS                     S32   2
@@ -407 +416 @@
 
 EXT.ANALYSIS.MAG.LIMITS METADATA
-    TYPE  DATA FILTER.ID MAG.LIMIT
-    gband DATA g         NAN            
-    rband DATA r         NAN             
-    iband DATA i         NAN             
-    zband DATA z         NAN             
-    yband DATA y         NAN             
-    wband DATA w         NAN             
-    other DATA any       NAN             
+    TYPE  DATA FILTER.ID MAG.LIMIT.PETRO    MAG.LIMIT.EXTFIT
+    gband DATA g         25                 22
+    rband DATA r         25                 22 
+    iband DATA i         25                 22
+    zband DATA z         25                 21
+    yband DATA y         25                 20
+    wband DATA w         25                 22
+    other DATA any       25                 22
 END
 
@@ -423 +432 @@
 RADIAL_APERTURES_SN_LIM             F32   0.0  # S/N limit for radial aperture calculation
 
+# zero point and exposure time to which fullForceSummary scales fluxes to
+PSPHOT.FULLFORCE.EXPTIME            F32   1.00 
+PSPHOT.FULLFORCE.ZERO_PT            F32   25.00
+
+
 GALAXY_SHAPES                       BOOL  F
+GALAXY_SHAPES_STYLE                 STR   R_MAJ_R_MIN
+
 GALAXY_SHAPES_FR_MAJOR_MIN          F32   0.5
-GALAXY_SHAPES_FR_MAJOR_MAX          F32   2.0
-GALAXY_SHAPES_FR_MAJOR_DEL          F32   0.1
+GALAXY_SHAPES_FR_MAJOR_MAX          F32   2.00
+GALAXY_SHAPES_FR_MAJOR_DEL          F32   0.05
 GALAXY_SHAPES_FR_MINOR_MIN          F32   0.5
-GALAXY_SHAPES_FR_MINOR_MAX          F32   2.0
-GALAXY_SHAPES_FR_MINOR_DEL          F32   0.1
-
+GALAXY_SHAPES_FR_MINOR_MAX          F32   2.00
+GALAXY_SHAPES_FR_MINOR_DEL          F32   0.05
+
+#GALAXY_SHAPES_FR_MAJOR_MIN          F32   0.85
+#GALAXY_SHAPES_FR_MAJOR_MAX          F32   1.15
+#GALAXY_SHAPES_FR_MAJOR_DEL          F32   0.05
+#GALAXY_SHAPES_FR_MINOR_MIN          F32   0.85
+#GALAXY_SHAPES_FR_MINOR_MAX          F32   1.15
+#GALAXY_SHAPES_FR_MINOR_DEL          F32   0.05
 
 # Extended source fit parameters
@@ -439 +461 @@
   EXTENDED_SOURCE_ANNULI              BOOL  TRUE
   EXT.NSIGMA.LIMIT.USE                BOOL  FALSE
+  EXT.FIT.MIN.GAL.LIMIT.USE           BOOL  TRUE  # limit extended source fits by galactic coordinates
 
   PSPHOT.STACK.MATCH.PSF.SOURCE       STR   AUTO # which inputs to convolve? (RAW, CNV, AUTO)
@@ -465 +488 @@
   EXT_MODEL                           STR   PS_MODEL_QGAUSS
   PEAKS_NMAX_TOTAL                    S32   0 # set this to limit the allowed number of peaks - Yields fault instead of avoid memory explosion
+
 
   SAVE.RESID                          BOOL  TRUE
@@ -683 +707 @@
     RADIAL_APERTURES        BOOL    TRUE
     EXTENDED_SOURCE_FITS    BOOL    TRUE    # this casues the xfit extension to be written out
-    EXTENDED_SOURCE_PETROSIAN BOOL    TRUE  # I want petrosian mags
+    EXTENDED_SOURCE_PETROSIAN BOOL  TRUE    # To measure petrosian magnitudes
+    EXTENDED_SOURCE_ANNULI  BOOL    TRUE
     SAVE.PSF                BOOL    FALSE
     SAVE.BACKMDL            BOOL    FALSE
-    # SAVE.RESID            BOOL    FALSE
-    # OUTPUT.FORMAT         STR     PS1_V3
-END
+    SAVE.RESID              BOOL    TRUE
+    OUTPUT.FORMAT           STR     PS1_SV2
+END

trunk/ppSub/src/ppSubMatchPSFs.c

@@ -156 +156 @@
     psLogMsg("ppSub", PS_LOG_INFO, "Input FWHM: %f\nReference FWHM: %f\n", inFWHM, refFWHM);
     if (!isfinite(inFWHM) || !isfinite(refFWHM)) {
-        psError(PPSUB_ERR_DATA, false, "Cannot determine FHWM for images, giving up.");
-        return false;
+        psErrorStackPrint(stderr, "Cannot determine FHWM for images, giving up.");
+        int error = psErrorCodeLast(); // Error code
+        ppSubDataQuality(data, error, PPSUB_FILES_ALL);
+        return true;
     }
 

trunk/psModules/src/camera/pmFPAfile.c

@@ -565 +565 @@
       case PM_FPA_FILE_CMF:
         return ("CMF");
+      case PM_FPA_FILE_CFF:
+        return ("CFF");
       case PM_FPA_FILE_WCS:
         return ("WCS");

trunk/psModules/src/objects/pmModelClass.c

@@ -66 +66 @@
 
 static pmModelClass *models = NULL;
+static psVector *modelClassLookupTable = NULL;  // translation between model types in header and here
 static int Nmodels = 0;
 
@@ -135 +136 @@
     models = NULL;
     Nmodels = 0;
+    psFree(modelClassLookupTable);
+    modelClassLookupTable = NULL;
     return;
 }
@@ -193 +196 @@
 }
 
+
+bool pmModelClassWriteHeader(psMetadata *header)
+{
+    psMetadataAddS32(header, PS_LIST_TAIL, "MTNUM", PS_META_REPLACE, "number of model types", Nmodels);
+    for (int i = 0; i < Nmodels; i++) {
+        char modelNameKey[16];
+        char modelValKey[16];
+        sprintf(modelNameKey, "MTNAM%02d", i);
+        sprintf(modelValKey,  "MTVAL%02d", i);
+        psMetadataAddStr(header, PS_LIST_TAIL, modelNameKey, PS_META_REPLACE, "", models[i].name);
+        psMetadataAddS32(header, PS_LIST_TAIL, modelValKey, PS_META_REPLACE, "", i);
+    }
+
+    return true;
+}
+
+bool pmModelClassReadHeader(psMetadata *header) {
+    psFree(modelClassLookupTable);
+
+    bool status;
+    int numHeaderModels = psMetadataLookupS32(&status, header, "MTNUM");
+    if (!status) {
+        return false;
+    }
+
+    psVector *inputTypes = psVectorAlloc(numHeaderModels, PS_TYPE_S32);
+    psVector *localTypes = psVectorAlloc(numHeaderModels, PS_TYPE_S32);
+    int max_val = -1;
+    for (int i = 0; i < numHeaderModels; i++) {
+        char modelNameKey[16];
+        char modelValKey[16];
+        sprintf(modelNameKey, "MTNAM%02d", i);
+        sprintf(modelValKey,  "MTVAL%02d", i);
+        psString thisName = psMetadataLookupStr(&status, header, modelNameKey);
+        int thisVal = psMetadataLookupS32(&status, header, modelValKey);
+        if (thisVal > max_val) {
+            max_val = thisVal;
+        }
+        inputTypes->data.S32[i] = thisVal;
+        localTypes->data.S32[i] = pmModelClassGetType(thisName);
+    }
+    if (max_val < 0) {
+        psFree(inputTypes);
+        psFree(localTypes);
+        return false;
+    }
+
+    modelClassLookupTable = psVectorAlloc(max_val + 1, PS_TYPE_S32);
+    psVectorInit(modelClassLookupTable, -1);
+
+    for (int i = 0; i < numHeaderModels; i++) {
+        int thisVal = inputTypes->data.S32[i];
+        int localVal = localTypes->data.S32[i];
+        modelClassLookupTable->data.S32[thisVal] = localVal;
+    }
+    psFree(inputTypes);
+    psFree(localTypes);
+
+    return true;
+}
+
+pmModelType pmModelClassGetLocalType(pmModelType inputType) {
+    pmModelType localType = -1;
+
+    if (modelClassLookupTable) {
+        if (inputType >= 0 && inputType < modelClassLookupTable->n) {
+            localType = modelClassLookupTable->data.S32[inputType];
+        }
+    } else {
+        // no lookup table defined
+        // for backwards compatability if inputType refers to a defined model, return it
+        if (inputType >= 0 && pmModelClassGetName(inputType)) {
+            localType = inputType;
+        }
+    }
+
+    return localType;
+}

trunk/psModules/src/objects/pmModelClass.h

@@ -76 +76 @@
 void pmModelClassSetLimits(pmModelLimitsType type);
 
+// write keywords to header definining the model type values used by this program
+bool pmModelClassWriteHeader(psMetadata *header);
+// create a lookup table for translating input model type values to local model type values
+bool pmModelClassReadHeader(psMetadata *header);
+// translate input model type value to local value
+pmModelType pmModelClassGetLocalType(pmModelType inputType);
 
 /// @}

trunk/psModules/src/objects/pmSourceIO.c

@@ -61 +61 @@
 static bool pmReadoutReadXFIT(pmFPAfile *file, pmReadout *readout, char * exttype, psMetadata *hduHeader, psString xfitname, psArray *sources, long *sourceIndex);
 static bool pmReadoutReadXRAD(pmFPAfile *file, pmReadout *readout, char * exttype, psMetadata *hduHeader, psString xfitname, psArray *sources, long *sourceIndex);
+static bool pmReadoutReadXGAL(pmFPAfile *file, pmReadout *readout, char * exttype, psMetadata *hduHeader, psString xfitname, psArray *sources, long *sourceIndex);
 
 // lookup the EXTNAME values used for table data and image header segments
@@ -374 +375 @@
         }                                                               \
         if (xgalname) {                                                 \
-            status &= pmSourcesWrite_##TYPE##_XGAL (file->fits, sources, xgalname, recipe); \
+            status &= pmSourcesWrite_##TYPE##_XGAL (file->fits, readout, sources, xgalname, recipe); \
         }                                                               \
     }
@@ -1036 +1037 @@
         bool XFIT_OUTPUT = psMetadataLookupBool(&status, recipe, "EXTENDED_SOURCE_FITS");
         bool XRAD_OUTPUT = psMetadataLookupBool(&status, recipe, "RADIAL_APERTURES");
-        bool XGAL_OUTPUT = false; // psMetadataLookupBool(&status, recipe, "GALAXY_SHAPES");
+        bool XGAL_OUTPUT = psMetadataLookupBool(&status, recipe, "GALAXY_SHAPES");
 
         if (!pmSourceIOextnames(&headname, &dataname, &deteffname, 
@@ -1122 +1123 @@
 
             long *sourceIndex = NULL;
-            if (XSRC_OUTPUT || XFIT_OUTPUT || XRAD_OUTPUT) {
+            if (XSRC_OUTPUT || XFIT_OUTPUT || XRAD_OUTPUT || XGAL_OUTPUT) {
+                // Build sourceIndex. Lookup table from source->seq to index in sources array.
+                // Consists of an array of length max(source->seq) + 1.
+
+                // find maximum sequence number
                 long seq_max = -1;
                 for (long i = sources->n -1; i >= 0; i--) {
@@ -1135 +1140 @@
                     }
                 }
+                // allocate and initialize the index
                 sourceIndex = psAlloc((seq_max + 1) * sizeof(long));
                 for (long i = 0; i < seq_max; i++) {
                     sourceIndex[i] = -1;
                 }
+                // populate the index
                 for (long i = 0; i < sources->n; i++) {
                     pmSource *source = sources->data[i];
@@ -1165 +1172 @@
                 psFree(xradname);
             }
+            if (XGAL_OUTPUT && xgalname) {
+                // a cmf file may have an XGAL extension, but it is not required
+                if (!pmReadoutReadXGAL(file, readout, exttype, hdu->header, xgalname, sources, sourceIndex)) {
+                    // do anything?
+                }
+                psFree(xgalname);
+            }
             psFree(sourceIndex);
 
@@ -1461 +1475 @@
     return status;
 }
+static bool pmReadoutReadXGAL(pmFPAfile *file, pmReadout *readout, char *exttype, psMetadata *hduHeader, psString xgalname, psArray *sources, long *sourceIndex) 
+{
+    if (!psFitsMoveExtNameClean (file->fits, xgalname)) {
+        psTrace ("pmFPAfile", 1, "cannot find xgal extension %s in %s, skipping", xgalname, file->filename);
+        return false;
+    }
+
+    psMetadata *tableHeader = psFitsReadHeader(NULL, file->fits); // The FITS header
+    if (!tableHeader) psAbort("cannot read table header");
+
+    char *xtension = psMetadataLookupStr (NULL, tableHeader, "XTENSION");
+    if (!xtension) psAbort("cannot read table type");
+    if (strcmp (xtension, "BINTABLE")) {
+        psFree(tableHeader);
+        psWarning ("no binary table in extension %s, skipping\n", xgalname);
+        return false;
+    }
+
+# define PM_SOURCES_READ_XGAL(NAME,TYPE)                                \
+    if (!strcmp (exttype, NAME)) {                                      \
+        status = pmSourcesRead_##TYPE##_XGAL(file->fits, readout, hduHeader, tableHeader, sources, sourceIndex); \
+    }                                                                   
+
+    bool status = false;
+    if (file->type == PM_FPA_FILE_CMF) {
+        PM_SOURCES_READ_XGAL("PS1_V1",    CMF_PS1_V1);
+        PM_SOURCES_READ_XGAL("PS1_V2",    CMF_PS1_V2);
+        PM_SOURCES_READ_XGAL("PS1_V3",    CMF_PS1_V3);
+        PM_SOURCES_READ_XGAL("PS1_V4",    CMF_PS1_V4);
+        PM_SOURCES_READ_XGAL("PS1_SV1",   CMF_PS1_SV1);
+        PM_SOURCES_READ_XGAL("PS1_SV2",   CMF_PS1_SV2);
+        PM_SOURCES_READ_XGAL("PS1_DV1",   CMF_PS1_DV1);
+        PM_SOURCES_READ_XGAL("PS1_DV2",   CMF_PS1_DV2);
+        PM_SOURCES_READ_XGAL("PS1_DV3",   CMF_PS1_DV3);
+    }
+    psFree(tableHeader);
+    return status;
+}

trunk/psModules/src/objects/pmSourceIO.h

@@ -21 +21 @@
   bool pmSourcesWrite_##TYPE##_XFIT(psFits *fits, pmReadout *readout, psArray *sources, psMetadata *imageHeader, char *extname); \
   bool pmSourcesWrite_##TYPE##_XRAD(psFits *fits, pmReadout *readout, psArray *sources, psMetadata *imageHeader, char *extname, psMetadata *recipe); \
-  bool pmSourcesWrite_##TYPE##_XGAL(psFits *fits, psArray *sources, char *extname, psMetadata *recipe); \
+  bool pmSourcesWrite_##TYPE##_XGAL(psFits *fits, pmReadout *readout, psArray *sources, char *extname, psMetadata *recipe); \
   psArray *pmSourcesRead_##TYPE (psFits *fits, psMetadata *header); \
   bool pmSourcesRead_##TYPE##_XSRC (psFits *fits, pmReadout *readout, psMetadata *header, psMetadata *tableHeader, psArray *sources, long *index); \
   bool pmSourcesRead_##TYPE##_XFIT (psFits *fits, pmReadout *readout, psMetadata *header, psMetadata *tableHeader, psArray *sources, long *index); \
   bool pmSourcesRead_##TYPE##_XRAD (psFits *fits, pmReadout *readout, psMetadata *header, psMetadata *tableHeader, psArray *sources, long *index);\
+  bool pmSourcesRead_##TYPE##_XGAL (psFits *fits, pmReadout *readout, psMetadata *header, psMetadata *tableHeader, psArray *sources, long *index);\
   
 // All of these functions need to use the same API, even if not all elements are used in a specific case

trunk/psModules/src/objects/pmSourceIO_CFF.c

@@ -65 +65 @@
     int modelType = pmModelClassGetType ("PS_MODEL_GAUSS");
 
+    // Read lookup table for model classes (if defined)
+    pmModelClassReadHeader(header);
+
     char *PSF_NAME = psMetadataLookupStr (&status, header, "PSFMODEL");
     if (PSF_NAME != NULL) {
@@ -111 +114 @@
         float theta      = psMetadataLookupF32 (&status, row, "THETA");
 
-        // XXX: we need to put a lookup table in the cff header to define the correspondence of the
-        // model type values in the cff with our models. (We want to use an interger for efficiency
-        // but the value for each model type is set on the organization of the the array in pmModelClass.c
-        // For now use the input values verbatim and trust the user that this is valid value
         int   galaxyModelType = psMetadataLookupS32(&status, row, "MODEL_TYPE");
+        if (status) {
+            galaxyModelType = pmModelClassGetLocalType(galaxyModelType);
+        } else {
+            galaxyModelType = -1;
+        }
         float Sindex     = psMetadataLookupF32 (&status, row, "INDEX"); // Should this be PAR_07 not sersic index
 
@@ -123 +127 @@
         source->type = PM_SOURCE_TYPE_STAR; // XXX this should be added to the flags
 
-        // XXX we can set this in general, but for a specific image, we need to weed out SATSTARS
+        // XXX we can set this in general, but for a specific image, we need to weed out SATSTARS and
+        // stars that are masked
         if (psfStar) {
             source->tmpFlags |= PM_SOURCE_TMPF_CANDIDATE_PSFSTAR;
@@ -180 +185 @@
         }
 
-        // XXX: should use < 0 as invalid galaxyModelType
-
-        if (fitGalaxy && galaxyModelType > 0) {
+        if (fitGalaxy && galaxyModelType >= 0) {
             source->modelFits = psArrayAllocEmpty (1);
             pmModel *model = pmModelAlloc(galaxyModelType);
@@ -235 +238 @@
     PS_ASSERT(mdok, false);
 
+    // write the definition of the model class type values to the header
+    psMetadata *outputHeader = psMetadataAlloc();
+    pmModelClassWriteHeader(outputHeader);
+
     psArray *table = psArrayAllocEmpty(sources->n);
 
@@ -246 +253 @@
         psS32 modelType = 0;
         bool fitGalaxy = false;
-        bool psfStar = false;
+        bool psfStar = (source->mode & PM_SOURCE_MODE_PSFSTAR) ? true : false;
         psF32 sersicIndex = 0;
-        if (source->modelFits == NULL) {
+        // For now only perform galaxy fits on extended objects
+        if (source->modelEXT == NULL) {
             pmModel *model = source->modelPSF;
             if (model == NULL) continue;
@@ -260 +268 @@
             yPos = model->params->data.F32[PM_PAR_YPOS];
             flux = source->psfFlux;
-            psfStar = (source->mode & PM_SOURCE_MODE_PSFSTAR) ? true : false;
             rMajor = 0;
             rMinor = 0;
@@ -329 +336 @@
 
         psArrayAdd(table, 100, row);
+        psFree(row);
     }
 
-    if (!psFitsWriteTable(fits, NULL, table, extname)) {
+    if (!psFitsWriteTable(fits, outputHeader, table, extname)) {
         psError(psErrorCodeLast(), false, "writing ext data %s\n", extname);
         psFree(table);
-        psFree(header);
+        psFree(outputHeader);
         return false;
     }
     psFree(table);
-    // psFree(header);
+    psFree(outputHeader);
 
     return true;

trunk/psModules/src/objects/pmSourceIO_CMF.c.in

@@ -708 +708 @@
             return false;
         }
-        // Find the source with this sequence number. 
-        // XXX: I am assuming that sources is sorted in order of seq
+        // Find the source with this sequence number using the sourceIndex. 
         long seq = psMetadataLookupU32 (&status, row, "IPP_IDET");
-        pmSource *source = NULL;
-#ifndef ASSUME_SORTED
-        long j = seq < sources->n ? seq : sources->n - 1;
-        for (; j >= 0; j--) {
-            source = sources->data[j];
-            if (source->seq == seq) {
-                break;
-            }
-        }
-#else
         long j = sourceIndex[seq];
         psAssert(j >= 0 && j < sources->n, "invalid sourceIndex");
-        source = sources->data[j];
-#endif
+        pmSource *source = sources->data[j];
         if (!source) {
             psError(PS_ERR_UNKNOWN, false, "Failed to find source for row %ld sequence number %ld\n", i, seq);
@@ -793 +781 @@
     // create a header to hold the output data
     psMetadata *outhead = psMetadataAlloc ();
+
+    pmModelClassWriteHeader(outhead);
 
     // write the links to the image header
@@ -1034 +1024 @@
         return false;
     }
+    // set up the lookup table to translate between input model types and output model types
+    // if not defined it is assumed that the tables are the same
+    pmModelClassReadHeader(tableHeader);
 
     for (long i = 0; i < numSources; i++) {
@@ -1042 +1035 @@
             return false;
         }
-        // Find the source with this sequence number. 
-        // XXX: I am assuming that sources is sorted in order of seq.
         long seq = psMetadataLookupU32 (&status, row, "IPP_IDET");
-        long j = seq < sources->n ? seq : sources->n - 1;
-        pmSource *source = NULL;
-        for (; j >= 0; j--) {
-            source = sources->data[j];
-            if (source->seq == seq) {
-                break;
-            }
-        }
+        long j = sourceIndex[seq];
+        psAssert(j >= 0 && j < sources->n, "invalid sourceIndex");
+        pmSource *source = sources->data[j];
         if (!source) {
             psError(PS_ERR_UNKNOWN, false, "Failed to find source for row %ld sequence number %ld\n", i, seq);
@@ -1095 +1081 @@
         // in the psf table.
         psS32 extModelType = psMetadataLookupS32(&status, row, "EXT_MODEL_TYPE");
-        if (!status) {
+        if (status) {
+            // translate between the type value in xfit and values used by this program
+            extModelType = pmModelClassGetLocalType(extModelType);
+        } else {
             // older cmfs don't have this column
             extModelType = -1;
-        }
+        } 
 
         psEllipseAxes axes;
@@ -1339 +1328 @@
             return false;
         }
-        // Find the source with this sequence number. 
-        // XXX: I am assuming that sources is sorted in order of seq.
         long seq = psMetadataLookupU32 (&status, row, "IPP_IDET");
-        long j = seq < sources->n ? seq : sources->n - 1;
-        pmSource *source = NULL;
-        for (; j >= 0; j--) {
-            source = sources->data[j];
-            if (source->seq == seq) {
-                break;
-            }
-        }
+        long j = sourceIndex[seq];
+        psAssert(j >= 0 && j < sources->n, "invalid sourceIndex");
+        pmSource *source = sources->data[j];
         if (!source) {
             psError(PS_ERR_UNKNOWN, false, "Failed to find source for row %ld sequence number %ld\n", i, seq);
@@ -1407 +1389 @@
 
 // XXX where should I record the number of columns??
-bool pmSourcesWrite_CMF_@CMFMODE@_XGAL (psFits *fits, psArray *sources, char *extname, psMetadata *recipe)
+bool pmSourcesWrite_CMF_@CMFMODE@_XGAL (psFits *fits, pmReadout *readout, psArray *sources, char *extname, psMetadata *recipe)
 {
     bool status = false;
@@ -1422 +1404 @@
     // write the links to the image header
     psMetadataAddStr (outhead, PS_LIST_TAIL, "EXTNAME", PS_META_REPLACE, "galaxy table extension", extname);
+
+    psMetadataAddStr (outhead, PS_LIST_TAIL, "HI", PS_META_REPLACE, "does this get through?", "THERE");
 
     // let's write these out in S/N order
@@ -1491 +1475 @@
 }
 
+bool pmSourcesRead_CMF_@CMFMODE@_XGAL(psFits *fits, pmReadout *readout, psMetadata *hduHeader, psMetadata *tableHeader, psArray *sources, long *sourceIndex)
+{
+    PS_ASSERT_PTR_NON_NULL(fits, false);
+    PS_ASSERT_PTR_NON_NULL(sources, false);
+
+    bool status;
+    long numSources = psFitsTableSize(fits); // Number of sources in table
+    if (numSources == 0) {
+        psError(psErrorCodeLast(), false, "XGAL Table contains no entries\n");
+        return false;
+    }
+
+    for (long i = 0; i < numSources; i++) {
+        psMetadata *row = psFitsReadTableRow(fits, i); // Table row
+        if (!row) {
+            psError(psErrorCodeLast(), false, "Unable to read row %ld of sources", i);
+            psFree(row);
+            return false;
+        }
+        // Find the source with this sequence number. 
+        // XXX: I am assuming that sources is sorted in order of seq
+        long seq = psMetadataLookupU32 (&status, row, "IPP_IDET");
+        long j = sourceIndex[seq];
+        psAssert(j >= 0 && j < sources->n, "invalid sourceIndex");
+
+        pmSource *source = sources->data[j];
+        if (!source) {
+            psError(PS_ERR_UNKNOWN, false, "Failed to find source for row %ld sequence number %ld\n", i, seq);
+            psFree(row);
+            return false;
+        }
+
+        psVector *Flux  = psMetadataLookupVector(&status, row, "GAL_FLUX");
+        psVector *dFlux = psMetadataLookupVector(&status, row, "GAL_FLUX_ERR");
+        psVector *chisq = psMetadataLookupVector(&status, row, "GAL_CHISQ");
+
+        if (Flux && Flux->n > 0) {
+            psFree(source->galaxyFits);
+            source->galaxyFits = pmSourceGalaxyFitsAlloc();
+            source->galaxyFits->nPix = psMetadataLookupF32(&status, row, "NPIX");
+
+            psFree(source->galaxyFits->Flux);
+            source->galaxyFits->Flux  = psMemIncrRefCounter(Flux);
+            psFree(source->galaxyFits->dFlux);
+            source->galaxyFits->dFlux = psMemIncrRefCounter(dFlux);
+            psFree(source->galaxyFits->chisq);
+            source->galaxyFits->chisq = psMemIncrRefCounter(chisq);
+        }
+
+        psFree(row);
+    }
+
+    return true;
+}

trunk/psModules/src/objects/pmSourceIO_PS1_CAL_0.c

@@ -714 +714 @@
 }
 
-bool pmSourcesWrite_PS1_CAL_0_XGAL (psFits *fits, psArray *sources, char *extname, psMetadata *recipe)
+bool pmSourcesWrite_PS1_CAL_0_XGAL (psFits *fits, pmReadout *readout, psArray *sources, char *extname, psMetadata *recipe)
 {
     return true;

trunk/psModules/src/objects/pmSourceIO_PS1_DEV_0.c

@@ -256 +256 @@
 }
 
-bool pmSourcesWrite_PS1_DEV_0_XGAL(psFits *fits, psArray *sources, char *extname, psMetadata *recipe)
-{
-    return true;
-}
+bool pmSourcesWrite_PS1_DEV_0_XGAL(psFits *fits, pmReadout *readout, psArray *sources, char *extname, psMetadata *recipe)
+{
+    return true;
+}

trunk/psModules/src/objects/pmSourceIO_PS1_DEV_1.c

@@ -596 +596 @@
 }
 
-bool pmSourcesWrite_PS1_DEV_1_XGAL(psFits *fits, psArray *sources, char *extname, psMetadata *recipe)
+bool pmSourcesWrite_PS1_DEV_1_XGAL(psFits *fits, pmReadout *readout, psArray *sources, char *extname, psMetadata *recipe)
 {
     return true;

trunk/psModules/src/objects/pmSourceIO_SMPDATA.c

@@ -226 +226 @@
 } 
 
-bool pmSourcesWrite_SMPDATA_XGAL(psFits *fits, psArray *sources, char *extname, psMetadata *recipe)
+bool pmSourcesWrite_SMPDATA_XGAL(psFits *fits, pmReadout *readout, psArray *sources, char *extname, psMetadata *recipe)
 {
     return true;

trunk/psastro/src/psastroDataSave.c

@@ -42 +42 @@
     }
 
-    // de-activate all files except PSASTRO.OUTPUT
+    // de-activate all files except PSASTRO.OUTPUT, PSASTRO.OUT.ASTROM, and PSPHOT.OUTPUT.CFF
     pmFPAfileActivate (config->files, false, NULL);
     pmFPAfileActivate (config->files, true, "PSASTRO.OUTPUT");
     pmFPAfileActivate (config->files, true, "PSASTRO.OUT.ASTROM");
+    pmFPAfileActivate (config->files, true, "PSPHOT.OUTPUT.CFF");
 
     pmFPAview *view = pmFPAviewAlloc (0);
@@ -103 +104 @@
     }
 
-    // activate all files except PSASTRO.OUTPUT
+    // activate all files except PSASTRO.OUTPUT, and PSPHOT.OUTPUT.CFF
     pmFPAfileActivate (config->files, true, NULL);
     pmFPAfileActivate (config->files, false, "PSASTRO.OUTPUT");
+    pmFPAfileActivate (config->files, false, "PSPHOT.OUTPUT.CFF");
 
     psFree (view);

trunk/psastro/src/psastroDefineFiles.c

@@ -100 +100 @@
     }
 
+    bool writeCff = psMetadataLookupBool (&status, recipe, "PSASTRO.SAVE.CFF");
+    if (writeCff) {
+        pmFPAfile *file = pmFPAfileDefineOutputFromFile  (config, input, "PSPHOT.OUTPUT.CFF");
+        if (!file) {
+            psError (PS_ERR_IO, false, "Can't find the output cff file definition");
+            return NULL;
+        }
+        if (file->type != PM_FPA_FILE_CFF) {
+            psError(PS_ERR_IO, true, "%s is not of type %s", "PSPHOT.OUTPUT.CFF", pmFPAfileStringFromType (PM_FPA_FILE_CFF));
+            return NULL;
+        }
+        file->save = true;
+    }
+
 
 # if (0)

trunk/psphot/src

@@ -26 +26 @@
 psphotFullForce
 psmakecff
+psphotFullForceSummary

@@ -26 +26 @@
 psphotFullForce
 psmakecff
+psphotFullForceSummary

trunk/psphot/src/Makefile.am

@@ -25 +25 @@
 libpsphot_la_LDFLAGS = $(PSPHOT_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
 
-bin_PROGRAMS = psphot psphotForced psphotFullForce psphotMinimal psphotMakePSF psphotStack psphotModelTest psmakecff
+bin_PROGRAMS = psphot psphotForced psphotFullForce psphotFullForceSummary psphotMinimal psphotMakePSF psphotStack psphotModelTest psmakecff
 # bin_PROGRAMS = psphotPetrosianStudy psphotTest psphotMomentsStudy 
 
@@ -39 +39 @@
 psphotFullForce_LDFLAGS = $(PSPHOT_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
 psphotFullForce_LDADD = libpsphot.la
+
+psphotFullForceSummary_CFLAGS = $(PSPHOT_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS)
+psphotFullForceSummary_LDFLAGS = $(PSPHOT_LIBS) $(PSMODULE_LIBS) $(PSLIB_LIBS)
+psphotFullForceSummary_LDADD = libpsphot.la
 
 psphotMinimal_CFLAGS = $(PSPHOT_CFLAGS) $(PSMODULE_CFLAGS) $(PSLIB_CFLAGS)
@@ -99 +103 @@
         psphotMosaicChip.c         \
         psphotCleanup.c
+
+# combine full force results from several inputs
+psphotFullForceSummary_SOURCES = \
+        psphotFullForceSummary.c \
+        psphotFullForceSummaryReadout.c
 
 # forced photometry of specified positions given a specified psf

trunk/psphot/src/psmakecff.c

@@ -2 +2 @@
 # include <config.h>
 # endif
+
+// psmakecff : A program to make read a cmf file and write a cff file
+// The real work is done in psModules.
 
 #include <stdio.h>
@@ -8 +11 @@
 #include "psphot.h"
 
+// For simplicilty, this program's (simple) functions are all contained in this file.
 static pmConfig* psmakecffArguments(int, char**);
 static bool psmakecffParseCamera(pmConfig *);
@@ -25 +29 @@
 //    psphotVersionPrint();
 
-    // load input data (config and images (signal, noise, mask)
     if (!psmakecffParseCamera (config)) {
         psErrorStackPrint(stderr, "Error setting up the camera\n");
@@ -31 +34 @@
     }
 
-    // call psphot for each readout
     if (!psmakecffImageLoop (config)) {
         psErrorStackPrint(stderr, "Error in the psphot image loop\n");
@@ -37 +39 @@
     }
 
-//    psLogMsg ("psphot", PS_LOG_WARN, "complete psphot run: %f sec\n", psTimerMark ("complete"));
-
-//    psErrorCode exit_status = psphotGetExitStatus();
-//    psphotCleanup (config);
     exit (0);
 }
@@ -73 +71 @@
     return config;
 }
+
 static bool psmakecffParseCamera(pmConfig *config) {
     bool status = false;
@@ -127 +126 @@
     psAssert (recipe, "missing recipe?");
 
-//    psImageMaskType maskTest = psMetadataLookupImageMask(&status, recipe, "MASK.PSPHOT");
 
     // for psphot, we force data to be read at the chip level
     while ((chip = pmFPAviewNextChip (view, input->fpa, 1)) != NULL) {
-        psLogMsg ("psmakecmf", 4, "Chip %d: %x %x\n", view->chip, chip->file_exists, chip->process);
+        psLogMsg ("psmakecff", 4, "Chip %d: %x %x\n", view->chip, chip->file_exists, chip->process);
         if (! chip->process || ! chip->file_exists) { continue; }
 
-#ifdef notdef
-        pmFPAfileActivate (config->files, false, NULL);
-        pmFPAfileActivate (config->files, true, "PSPHOT.LOAD");
-#endif
         if (!pmFPAfileIOChecks (config, view, PM_FPA_BEFORE)) ESCAPE ("failed input for Chip in psmakecff.");
 
         // there is now only a single chip (multiple readouts?). loop over it and process
         while ((cell = pmFPAviewNextCell (view, input->fpa, 1)) != NULL) {
-            psLogMsg ("psmakecmf", 5, "Cell %d: %x %x\n", view->cell, cell->file_exists, cell->process);
+            psLogMsg ("psmakecff", 5, "Cell %d: %x %x\n", view->cell, cell->file_exists, cell->process);
 
             // process each of the readouts
             while ((readout = pmFPAviewNextReadout (view, input->fpa, 1)) != NULL) {
-                psLogMsg ("psmakecmf", 6, "Cell %d: %x %x\n", view->cell, cell->file_exists, cell->process);
+                psLogMsg ("psmakecff", 6, "Cell %d: %x %x\n", view->cell, cell->file_exists, cell->process);
                 if (! readout->data_exists) { continue; }
 
                 pmHDU *hdu = pmHDUGetHighest(input->fpa, chip, cell);
                 if (hdu && hdu != lastHDU) {
+                    // XXX: probably should do this
                     // psphotVersionHeaderFull(hdu->header);
                     lastHDU = hdu;
@@ -157 +152 @@
 
         }
-        // Defer output and closing of files until we've (possibly) done the NFrames analysis below
-//        pmFPAfileActivate (config->files, false, NULL);
-
         if (!pmFPAfileIOChecks (config, view, PM_FPA_AFTER)) ESCAPE ("failed pmFPAfileIOChecks for Chip in psmakecff.");
     }
     if (!pmFPAfileIOChecks (config, view, PM_FPA_AFTER)) ESCAPE ("failed pmFPAfileIOChecks for FPA in psphot.");
 
-    // fail if we failed to handle an error
+    // fail if we encountered an unhandled error
     if (psErrorCodeLast() != PS_ERR_NONE) psAbort ("failed to handle an error!");
 

trunk/psphot/src/psphot.h

@@ -361 +361 @@
 bool psphotModelTestReadout(pmConfig *config, const pmFPAview *view, const char *filerule);
 
+bool psphotFullForceSummaryReadout (pmConfig * config, const pmFPAview *view);
+
 int psphotFileruleCount(const pmConfig *config, const char *filerule);
 bool psphotFileruleCountSet(const pmConfig *config, const char *filerule, int num);

trunk/psphot/src/psphotChooseAnalysisOptions.c

@@ -60 +60 @@
 }
 
-bool GetGalacticCoords (psSphere *ptGal, psSphere *ptSky, psSphereRot *toGal, pmChip *chip, float xPos, float yPos);
+static bool GetGalacticCoords (psSphere *ptGal, psSphere *ptSky, psSphereRot *toGal, pmChip *chip, float xPos, float yPos);
 
 // this function use an internal flag to mark sources which have already been measured
@@ -121 +121 @@
     psSphereRot *toGal = NULL;
     float GAL_LIMIT = 0.0;
+    float GAL_LIMIT_BULGE = 0.0;
+    float GAL_LIMIT_SIGMA2 = 0.0;
     bool useGAL_LIMIT = psMetadataLookupBool(&status, recipe, "EXT.FIT.MIN.GAL.LIMIT.USE");
     assert (status);
     if (useGAL_LIMIT) {
+        toGal = psSphereRotICRSToGalactic();
         GAL_LIMIT = psMetadataLookupF32(&status, recipe, "EXT.FIT.MIN.GAL.LIMIT");
-        toGal = psSphereRotICRSToGalactic();
+        assert (status);
+        GAL_LIMIT = DEG_TO_RAD(GAL_LIMIT);
+        GAL_LIMIT_BULGE = psMetadataLookupF32(&status, recipe, "EXT.FIT.MIN.GAL.LIMIT.BULGE");
+        assert (status);
+        GAL_LIMIT_BULGE = DEG_TO_RAD(GAL_LIMIT_BULGE);
+        float GAL_LIMIT_SIGMA = psMetadataLookupF32(&status, recipe, "EXT.FIT.MIN.GAL.LIMIT.BULGE.SIGMA");
+        assert (status);
+        assert(GAL_LIMIT_SIGMA > 0);
+        GAL_LIMIT_SIGMA = DEG_TO_RAD(GAL_LIMIT_SIGMA);
+        GAL_LIMIT_SIGMA2 = GAL_LIMIT_SIGMA * GAL_LIMIT_SIGMA;
     }
 
@@ -135 +147 @@
     }
 
+    float petroFluxLim = NAN;
     float extFitFluxLim = NAN;
 
@@ -140 +153 @@
         float extFitMagLimDefault = NAN;
         float extFitMagLim = NAN;
+        float petroMagLimDefault = NAN;
+        float petroMagLim = NAN;
 
         // match to the given filter
@@ -156 +171 @@
             // find a matching filter or default to 'any'
             if (!strcasecmp (thisFilter, "any")) {
-                psString extFitMagLimStr = psMetadataLookupStr (&status, item->data.md, "MAG.LIMIT");
-                psAssert(extFitMagLimStr, "missing MAG.LIMIT");
+                psString petroMagLimStr = psMetadataLookupStr (&status, item->data.md, "MAG.LIMIT.PETRO");
+                psAssert(petroMagLimStr, "missing MAG.LIMIT.PETRO");
+                petroMagLimDefault = atof (petroMagLimStr);
+
+                psString extFitMagLimStr = psMetadataLookupStr (&status, item->data.md, "MAG.LIMIT.EXTFIT");
+                psAssert(extFitMagLimStr, "missing MAG.LIMIT.EXTFIT");
                 extFitMagLimDefault = atof (extFitMagLimStr);
             }
 
-            // find a matching filter or default to 'any'
             if (!strcasecmp (thisFilter, filterID)) {
-                psString extFitMagLimStr = psMetadataLookupStr (&status, item->data.md, "MAG.LIMIT");
-                psAssert(extFitMagLimStr, "missing MAG.LIMIT");
+                psString petroMagLimStr = psMetadataLookupStr (&status, item->data.md, "MAG.LIMIT.PETRO");
+                psAssert(petroMagLimStr, "missing MAG.LIMIT.PETRO");
+                petroMagLim = atof (petroMagLimStr);
+
+                psString extFitMagLimStr = psMetadataLookupStr (&status, item->data.md, "MAG.LIMIT.EXTFIT");
+                psAssert(extFitMagLimStr, "missing MAG.LIMIT.EXTFIT");
                 extFitMagLim = atof (extFitMagLimStr);
                 break;
             }
         }
+        psFree(iter);
+        if (!isfinite (petroMagLim)) petroMagLim = petroMagLimDefault;
         if (!isfinite (extFitMagLim)) extFitMagLim = extFitMagLimDefault;
 
@@ -185 +209 @@
         psAssert (status, "missing FPA.ZP?");
 
+        petroFluxLim = exptime * pow (10.0, 0.4*(zeropt - petroMagLim));
         extFitFluxLim = exptime * pow (10.0, 0.4*(zeropt - extFitMagLim));
     }
@@ -225 +250 @@
             psSphere ptGal, ptSky;
             GetGalacticCoords (&ptGal, &ptSky, toGal, chip, source->peak->xf, source->peak->yf);
-            if (fabs(ptGal.d) < GAL_LIMIT) continue;
+            float b = ptGal.r;
+            float limit = GAL_LIMIT + GAL_LIMIT_BULGE * exp(-0.5*(b*b/GAL_LIMIT_SIGMA2));
+            if (fabs(ptGal.d) < limit) continue;
             // include an exception for low density skycells below the limit?
         }
 
         // for petro and extFit, we will either use the mag limits or the S/N
+        if (doPetrosian) {
+            if (isfinite(petroFluxLim)) {
+                if (source->moments->KronFlux > petroFluxLim) {
+                    source->tmpFlags |= PM_SOURCE_TMPF_PETRO;
+                }
+            } else if (source->moments->KronFlux > SN_LIM * source->moments->KronFluxErr) {
+                source->tmpFlags |= PM_SOURCE_TMPF_PETRO;
+            }
+        }
         if (isfinite(extFitFluxLim)) {
             if (source->moments->KronFlux > extFitFluxLim) {
                 source->tmpFlags |= PM_SOURCE_TMPF_EXT_FIT;
-                if (doPetrosian) source->tmpFlags |= PM_SOURCE_TMPF_PETRO;
             }
         } else {
             if (source->moments->KronFlux > SN_LIM * source->moments->KronFluxErr) {
                 source->tmpFlags |= PM_SOURCE_TMPF_EXT_FIT;
-                if (doPetrosian) source->tmpFlags |= PM_SOURCE_TMPF_PETRO;
             }
         }
@@ -244 +278 @@
 
     psLogMsg ("psphot.options", PS_LOG_WARN, "choose analysis options for %ld sources: %f sec\n", sources->n, psTimerMark ("psphot.options"));
+
+    psFree(toGal);
 
     return true;
@@ -287 +323 @@
     psSphereRot *toGal = NULL;
     float GAL_LIMIT = 0.0;
+    float GAL_LIMIT_BULGE = 0.0;
+    float GAL_LIMIT_SIGMA2 = 0.0;
     bool useGAL_LIMIT = psMetadataLookupBool(&status, recipe, "EXT.FIT.MIN.GAL.LIMIT.USE");
     assert (status);
     if (useGAL_LIMIT) {
+        toGal = psSphereRotICRSToGalactic();
         GAL_LIMIT = psMetadataLookupF32(&status, recipe, "EXT.FIT.MIN.GAL.LIMIT");
-        toGal = psSphereRotICRSToGalactic();
+        assert (status);
+        GAL_LIMIT = DEG_TO_RAD(GAL_LIMIT);
+        GAL_LIMIT_BULGE = psMetadataLookupF32(&status, recipe, "EXT.FIT.MIN.GAL.LIMIT.BULGE");
+        assert (status);
+        GAL_LIMIT_BULGE = DEG_TO_RAD(GAL_LIMIT_BULGE);
+        float GAL_LIMIT_SIGMA = psMetadataLookupF32(&status, recipe, "EXT.FIT.MIN.GAL.LIMIT.BULGE.SIGMA");
+        assert (status);
+        assert(GAL_LIMIT_SIGMA > 0);
+        GAL_LIMIT_SIGMA = DEG_TO_RAD(GAL_LIMIT_SIGMA);
+        GAL_LIMIT_SIGMA2 = GAL_LIMIT_SIGMA * GAL_LIMIT_SIGMA;
     }
 
@@ -365 +413 @@
     // find extFitFluxLim->data.F32[i] for i == image number
     psVector *extFitFluxLim = NULL;
+    psVector *petroFluxLim = NULL;
     if (magLimits) {
         extFitFluxLim = psVectorAlloc (inputFilters->n, PS_TYPE_F32);
+        petroFluxLim = psVectorAlloc (inputFilters->n, PS_TYPE_F32);
         psVectorInit (extFitFluxLim, NAN);
+        psVectorInit (petroFluxLim, NAN);
 
         float extFitMagLimDefault = NAN;
+        float petroMagLimDefault = NAN;
 
         // match mag limits (flux limits) to the filters
@@ -383 +435 @@
             // save the default value to assign to unset filters
             if (!strcasecmp (thisFilter, "any")) {
-                psString magString = psMetadataLookupStr (&status, item->data.md, "MAG.LIMIT");
-                psAssert(magString, "missing MAG.LIMIT");
+                psString magString = psMetadataLookupStr (&status, item->data.md, "MAG.LIMIT.EXTFIT");
+                psAssert(magString, "missing MAG.LIMIT.EXTFIT");
                 extFitMagLimDefault = atof (magString);
+
+                magString = psMetadataLookupStr (&status, item->data.md, "MAG.LIMIT.PETRO");
+                psAssert(magString, "missing MAG.LIMIT.PETRO");
+                petroMagLimDefault = atof (magString);
                 continue;
             }
@@ -393 +449 @@
                 if (i == chisqNum) continue;
                 if (!strcasecmp (thisFilter, inputFilters->data[i])) {
-                    psString magString = psMetadataLookupStr (&status, item->data.md, "MAG.LIMIT");
-                    psAssert(magString, "missing MAG.LIMIT");
+                    psString magString = psMetadataLookupStr (&status, item->data.md, "MAG.LIMIT.PETRO");
+                    psAssert(magString, "missing MAG.LIMIT.PETRO");
                     float magvalue = atof (magString);
+                    petroFluxLim->data.F32[i] = exptime->data.F32[i] * pow (10.0, 0.4*(zeropt->data.F32[i] - magvalue));
+
+                    magString = psMetadataLookupStr (&status, item->data.md, "MAG.LIMIT.EXTFIT");
+                    psAssert(magString, "missing MAG.LIMIT.EXTFIT");
+                    magvalue = atof (magString);
                     extFitFluxLim->data.F32[i] = exptime->data.F32[i] * pow (10.0, 0.4*(zeropt->data.F32[i] - magvalue));
                     break;
@@ -401 +462 @@
             }
         }
+        psFree(iter);
 
         for (int i = 0; i < num; i++) {
             if (i == chisqNum) continue;
-            if (isfinite(extFitFluxLim->data.F32[i])) continue;
-            extFitFluxLim->data.F32[i] = exptime->data.F32[i] * pow (10.0, 0.4*(zeropt->data.F32[i] - extFitMagLimDefault));
+            if (!isfinite(petroFluxLim->data.F32[i])) {
+                petroFluxLim->data.F32[i] = exptime->data.F32[i] * pow (10.0, 0.4*(zeropt->data.F32[i] - petroMagLimDefault));
+            }
+            if (!isfinite(extFitFluxLim->data.F32[i])) {
+                extFitFluxLim->data.F32[i] = exptime->data.F32[i] * pow (10.0, 0.4*(zeropt->data.F32[i] - extFitMagLimDefault));
+            }
         }
     }
@@ -420 +486 @@
         bool doObjectRadial = false;
         bool doObjectExtFit = false;
-        for (int j = 0; !doObjectExtFit && !doObjectRadial && (j < object->sources->n); j++) {
+        bool doObjectPetrosian = false;
+        for (int j = 0; !doObjectExtFit && !doObjectRadial && !doObjectPetrosian && (j < object->sources->n); j++) {
 
             pmSource *source = object->sources->data[j];
@@ -451 +518 @@
             if (extFitAll) {
                 doObjectExtFit = true;
+                doObjectPetrosian = doPetrosian;
                 continue;
             }
@@ -465 +533 @@
                 psSphere ptGal, ptSky;
                 GetGalacticCoords (&ptGal, &ptSky, toGal, chips->data[imageID], source->peak->xf, source->peak->yf);
-                if (fabs(ptGal.d) < GAL_LIMIT) continue;
-                // include an exception for low density skycells below the limit?
-            }
-
-            float fluxLim = extFitFluxLim ? extFitFluxLim->data.F32[imageID] : NAN;
+                float b = ptGal.r;
+                float limit = GAL_LIMIT + GAL_LIMIT_BULGE * exp(-0.5*(b*b/GAL_LIMIT_SIGMA2));
+                if (fabs(ptGal.d) < limit) continue;
+            }
+
+            float fluxLim = NAN;
+            // for petro and extFit, we will either use the mag limits or the S/N
+            if (doPetrosian) {
+                fluxLim = petroFluxLim ? petroFluxLim->data.F32[imageID] : NAN;
+                if (isfinite(fluxLim)) {
+                    if (source->moments->KronFlux > extFitFluxLim->data.F32[imageID]) {
+                        doObjectPetrosian = true;
+                    }
+                } else {
+                    if (source->moments->KronFlux > SN_LIM * source->moments->KronFluxErr) {
+                        doObjectPetrosian = true;
+                    }
+                }
+            }
+            fluxLim = extFitFluxLim ? extFitFluxLim->data.F32[imageID] : NAN;
             // for petro and extFit, we will either use the mag limits or the S/N
             if (isfinite(fluxLim)) {
@@ -503 +586 @@
             if (source->modelPSF == NULL) continue;
                 
-            // Do the fits if the recipe requests we do extended source fits to everything
             if (doObjectExtFit) {
                 source->tmpFlags |= PM_SOURCE_TMPF_EXT_FIT;
-                if (doPetrosian) source->tmpFlags |= PM_SOURCE_TMPF_PETRO;
-            }
+            }
+            if (doObjectPetrosian) {
+                source->tmpFlags |= PM_SOURCE_TMPF_PETRO;
+            }
 
             // Do the fits if the recipe requests we do extended source fits to everything
@@ -522 +606 @@
     psLogMsg ("psphot.options", PS_LOG_WARN, "choose analysis options for %ld objects: %f sec\n", objects->n, psTimerMark ("psphot.options"));
 
+    psFree(toGal);
+
     return true;
 }
 
-bool GetGalacticCoords (psSphere *ptGal, psSphere *ptSky, psSphereRot *toGal, pmChip *chip, float xPos, float yPos) {
+static bool GetGalacticCoords (psSphere *ptGal, psSphere *ptSky, psSphereRot *toGal, pmChip *chip, float xPos, float yPos) {
 
     pmFPA *fpa = chip->parent;
@@ -543 +629 @@
     psSphereRotApply (ptGal, toGal, ptSky);
 
+    // psSphereRotApply insures that 0 < r < 2PI. We want -PI < b <= PI
+    if (ptGal->r > M_PI) {
+        ptGal->r -= 2.0 * M_PI;
+    }
+
     return true;
 

trunk/psphot/src/psphotGalaxyShape.c

@@ -283 +283 @@
     }
 
+#ifdef SAVE_BEST_MODEL
+    // Save model with smallest chisq
     if (isfinite(chisqBest)) {
+#else 
+    // Save model with nominal parameters
+    {
+        fRmajorBest = 1;
+        fRminorBest = 1;
+#endif
         // now save the best fitting model as the source's extended model
         psEllipseAxes testAxes = guessAxes;
@@ -293 +301 @@
         psphotGalaxyShapeSource (pcm, source, maskVal, psfSize, false);
 
+        // Replace modelEXT with this model
+        // XXX: only do this if the model is good
         psFree (source->modelEXT);
 
@@ -299 +309 @@
         source->mode |= PM_SOURCE_MODE_EXTMODEL;
         source->mode |= PM_SOURCE_MODE_NONLINEAR_FIT;
-
-        // adjust the window so the subtraction covers the faint wings
-        // psphotSetRadiusMoments(&fitRadius, &windowRadius, readout, source, markVal);
 
         // cache the model flux

@@ -47 +47 @@
 sctool
 releasetool
+fftool

@@ -26 +26 @@
 psphotFullForce
 psmakecff
+psphotFullForceSummary