Index: branches/eam_branches/ipp-20140904/Nebulous/bin/neb-rm =================================================================== --- branches/eam_branches/ipp-20140904/Nebulous/bin/neb-rm (revision 37803) +++ branches/eam_branches/ipp-20140904/Nebulous/bin/neb-rm (revision 37804) @@ -44,5 +44,7 @@ foreach my $key (@ARGV) { unless (defined $move) { - $neb->delete($key, $force,$invalid) or warn "failed to delete key: $key"; + eval { + $neb->delete($key, $force,$invalid) or warn "failed to delete key: $key"; + } } else { $neb->move($key, "${key}.rm") or warn "failed to move key: $key"; Index: branches/eam_branches/ipp-20140904/extsrc/gpcsw/gpcsrc/fits/burntool/persistio.c =================================================================== --- branches/eam_branches/ipp-20140904/extsrc/gpcsw/gpcsrc/fits/burntool/persistio.c (revision 37803) +++ branches/eam_branches/ipp-20140904/extsrc/gpcsw/gpcsrc/fits/burntool/persistio.c (revision 37804) @@ -115,4 +115,7 @@ #define DIFFERENT_STREAK 20 /* Proximity for union versus intersection */ +# define myMAXSIZE 20000 + + /****************************************************************/ /* persist_merge(): Disentangle overlapping persistence streaks */ @@ -123,5 +126,5 @@ int lapmax; OBJBOX *box; - int xusage[MAXSIZE], yctr[MAXSIZE], boxid[MAXSIZE]; + int xusage[myMAXSIZE], yctr[myMAXSIZE], boxid[myMAXSIZE]; double zk, zj; @@ -131,15 +134,15 @@ /* Assess the clustering of all the streaks */ - for(i=0; i= MAXSIZE-1) continue; + if(box[k].exfit >= myMAXSIZE-1) continue; for(i=box[k].sxfit; i<=box[k].exfit; i++) xusage[i] += 1; } /* Identify clusters */ - for(xs=0; xsreduction($reduction, 'ADDSTAR'); # XXX This is now not used: do we still need it? @@ -132,4 +141,11 @@ print "using $fpaObjects for $stage\n"; } +if ($stage =~/diff/) { + print "using $fpaObjects for $stage\n"; +} +if ($stage =~/fullforce/) { + print "using $fpaObjects for $stage\n"; +} + @@ -177,5 +193,5 @@ } -if (($stage =~/stack/ || $stage =~ /skycal/)) { +if (($stage =~/stack/ || $stage =~ /skycal/ || $stage =~/diff/ || $stage=~/fullforce/)) { $fpaObjects =~ s/smf$/cmf/; $fpaObjectsAlt =~ s/smf$/cmf/; @@ -305,4 +321,10 @@ $command .= " -quick-airmass"; } #careful here - this matches staticsky and staticsky_multi + if ($stage =~ /diff/) { + $command .= " -accept-astrom"; + } + if ($stage =~ /fullforce/) { + $command .= " -accept-astrom -xrad"; + } my $mjd_addstar_start = DateTime->now->mjd; # MJD of starting script Index: branches/eam_branches/ipp-20140904/ippScripts/scripts/ipp_cleanup.pl =================================================================== --- branches/eam_branches/ipp-20140904/ippScripts/scripts/ipp_cleanup.pl (revision 37803) +++ branches/eam_branches/ipp-20140904/ippScripts/scripts/ipp_cleanup.pl (revision 37804) @@ -588,4 +588,5 @@ addFilename(\@files, "PSWARP.OUTPUT.SOURCES", $path_base, $skycell_id); addFilename(\@files, "SKYCELL.TEMPLATE", $path_base, $skycell_id ); + addFilename(\@files, "PSPHOT.RESID",$path_base,"fpa"); if ($mode eq "goto_purged") { # additional files to remove for 'purge' mode Index: branches/eam_branches/ipp-20140904/ippScripts/scripts/lap_science.pl =================================================================== --- branches/eam_branches/ipp-20140904/ippScripts/scripts/lap_science.pl (revision 37803) +++ branches/eam_branches/ipp-20140904/ippScripts/scripts/lap_science.pl (revision 37804) @@ -622,4 +622,8 @@ $complete_fstacks++; } + elsif (($stack->{final_state})&& + ($stack->{final_state} eq 'wait')) { + $complete_fstacks++; # Not the best either, but we've noted an issue, and tried to move on. LAP should respect this. + } elsif (($stack->{final_state})&&($stack->{final_state} eq 'new')&& ($stack->{final_fault} >= 4)&&($stack->{final_fault} != 32767)) { @@ -757,5 +761,5 @@ } elsif ((defined($exposure->{warpRun_state}))&&($exposure->{warpRun_state} eq 'full')&& - ($exposure->{warp_bad_quality} / $exposure->{warp_component_count} > 0.2)) { + ($exposure->{warp_bad_quality} / $exposure->{warp_component_count} > 0.3)) { printf("QUALITY: $exposure->{exp_id} has bad warp quality: %d / %d\n", $exposure->{warp_bad_quality} , $exposure->{warp_component_count}); Index: branches/eam_branches/ipp-20140904/ippScripts/scripts/permcheck.pl =================================================================== --- branches/eam_branches/ipp-20140904/ippScripts/scripts/permcheck.pl (revision 37803) +++ branches/eam_branches/ipp-20140904/ippScripts/scripts/permcheck.pl (revision 37804) @@ -21,4 +21,5 @@ #my $regtool = can_run('regtool') or (warn "Can't find regtool" and $missing_tools = 1); my $camtool = can_run('camtool') or (warn "Can't find camtool" and $missing_tools = 1); +my $warptool = can_run('warptool') or (warn "Can't find warptool" and $missing_tools = 1); my $stacktool = can_run('stacktool') or (warn "Can't find stacktool" and $missing_tools = 1); my $staticskytool = can_run('staticskytool') or (warn "Can't find staticskytool" and $missing_tools = 1); @@ -41,5 +42,5 @@ pod2usage( -msg => "Required options: --dbname", -exitval => 2 ) unless $dbname; -pod2usage( -msg => "Required options: --stage [camera|stack|skycal]", -exitval => 2 ) +pod2usage( -msg => "Required options: --stage [camera|warp|stack|skycal]", -exitval => 2 ) unless $stage; pod2usage( -msg => "missing key", exitval => 2 ) @@ -58,4 +59,15 @@ 'ippb06' => 1 ); + +if ($stage eq 'warp') { # These things are destined to have their only copy on the stsciXX nodes. + %backup_hosts = ('stsci00' => 1, 'stsci01' => 1, 'stsci02' => 1, + 'stsci03' => 1, 'stsci04' => 1, 'stsci05' => 1, 'stsci06' => 1, + 'stsci07' => 1, 'stsci08' => 1, 'stsci09' => 1, + 'stsci10' => 1, 'stsci11' => 1, 'stsci12' => 1, + 'stsci13' => 1, 'stsci14' => 1, 'stsci15' => 1, 'stsci16' => 1, + 'stsci17' => 1, 'stsci18' => 1, 'stsci19' => 1); + %backup_destinations = %backup_hosts; +} + my $backup_Nvols = 3; @@ -113,6 +125,18 @@ &my_die("Unable to parse metadata from stagetool", $stage_id); } +elsif ($stage eq 'warp') { + my $cmd = "$warptool -warped -warp_id $stage_id -dbname $dbname"; + my ($success, $error_code, $full_buf, $stdout_buf, $stderr_buf ) = + run(command => $cmd, verbose => 0); + unless ($success) { + $error_code = (($error_code >> 8) or $PS_EXIT_PROG_ERROR); + &my_die("Unable to perform stagetool: $error_code", $stage_id); + } + $files = $mdcParser->parse_list(join "", @$stdout_buf) or + &my_die("Unable to parse metadata from stagetool", $stage_id); +} my %components = ('camera' => ['PSASTRO.OUTPUT','PSASTRO.OUTPUT.MASK'], + 'warp' => ['PSWARP.OUTPUT','PSWARP.OUTPUT.MASK','PSWARP.OUTPUT.VARIANCE'], 'stack' => ['PPSTACK.UNCONV.COMP','PPSTACK.UNCONV.MASK.COMP','PPSTACK.UNCONV.VARIANCE.COMP', 'PPSTACK.UNCONV.EXP','PPSTACK.UNCONV.EXPNUM','PPSTACK.UNCONV.EXPWT.COMP'], @@ -124,4 +148,7 @@ my $data_state = $entry->{state}; my $hostname = $entry->{hostname}; + my $quality = $entry->{quality}; + if ($quality != 0) { next; } + print "# $path_base $data_state $hostname $quality\n"; foreach my $product (@{ $components{$stage} }) { my @keys = (); @@ -148,5 +175,9 @@ my $stat = $neb->stat($key); die "nebulous key: $key not found" unless $stat; - my $instances = $neb->find_instances($key, 'any'); + my $instances; + eval { + $instances = $neb->find_instances($key, 'any'); + }; + unless (defined($instances)) { print "## skipping due to zero instances\n"; next; } die "no instances found" unless $instances; @@ -256,5 +287,5 @@ $quality[0] = 1; print "cp $good_file $diskfiles[0]\n"; - system("cp $good_file $diskfiles[0]"); + vsystem("cp $good_file $diskfiles[0]"); } # Begin my best validation thought @@ -275,5 +306,5 @@ if ($md5sums[$i] eq $md5sum) { # Found it. print "cp $diskfiles[$i] $diskfiles[0]\n"; - system("cp $diskfiles[$i] $diskfiles[0]"); + vsystem("cp $diskfiles[$i] $diskfiles[0]"); # Begin my best validation thought { @@ -291,8 +322,10 @@ if ($Ngood == 1) { # We have only one version if ($quality_mask & 1) { # ANd it's on a backup volume - print "neb-replicate $key\n"; - if ($do_ops) { - $neb->replicate($key) or die "failed to replicate the single valid copy"; - if ($@) { die $@; } + if ($user_copies > 1) { # And we want more than one copy. + print "neb-replicate $key\n"; + if ($do_ops) { + $neb->replicate($key) or die "failed to replicate the single valid copy"; + if ($@) { die $@; } + } } } @@ -305,5 +338,5 @@ # Begin my best validation thought - system("sync") == 0 or die "Couldn't sync?"; + vsystem("sync") == 0 or die "Couldn't sync?"; my $uris = $neb->find_instances($key,$rep_vol); @$uris = map {URI->new($_)->file if $_} @$uris; @@ -325,5 +358,5 @@ # Begin my best validation thought - system("sync") == 0 or die "Couldn't sync?"; + vsystem("sync") == 0 or die "Couldn't sync?"; my $uris = $neb->find_instances($key,$rep_vol); @$uris = map {URI->new($_)->file if $_} @$uris; @@ -343,23 +376,36 @@ } } - for (my $i = 0; $i <= $#diskfiles; $i++) { + } +#XXX NEW DEBUG THIS: This should iterate over diskfiles that are not marked with a quality = 2, and cull them. +# For the N>=2, we set that for the primary, and +# for both cases, we do not have a diskfile entry for the newly replicated copy. +# Therefore, this should cull down to the correct number. + for (my $i = 0; $i <= $#diskfiles; $i++) { # print "$existance[$i] $quality[$i] $md5sums[$i] $md5sum $diskhosts[$i] $diskvols[$i]\n"; - if ($quality[$i] != 1) { - if ($existance[$i] == 0) { # This disk file doesn't exist. - system("touch $diskfiles[$i]"); - } - my $cull_vol = $diskhosts[$i] . "." . $diskvols[$i]; - print "neb-cull --volume $cull_vol $key\n"; - if ($do_ops) { - # The tilde here is to force hard volumes. Don't touch it. - $neb->cull($key,"~${cull_vol}",2) or die "failed to cull a superfluous instance"; - if ($@) { die "$@"; } - } - } - } - } + if ($quality[$i] != 1) { + if ($existance[$i] == 0) { # This disk file doesn't exist. + vsystem("touch $diskfiles[$i]"); + } + my $cull_vol = $diskhosts[$i] . "." . $diskvols[$i]; + print "neb-cull --volume $cull_vol $key\n"; + if ($do_ops) { + # The tilde here is to force hard volumes. Don't touch it. + $neb->cull($key,"~${cull_vol}",2) or die "failed to cull a superfluous instance"; + if ($@) { die "$@"; } + } + } + } + # } } # end keys for this product } # end product for this entry } # end entry for this id + +sub vsystem { + my $cmd = shift; + print "$cmd\n"; + if ($do_ops) { + system($cmd); + } +} sub local_md5sum { @@ -379,4 +425,5 @@ + sub parse_volume { my $filename = shift(@_); Index: branches/eam_branches/ipp-20140904/ippScripts/scripts/publish_file.pl =================================================================== --- branches/eam_branches/ipp-20140904/ippScripts/scripts/publish_file.pl (revision 37803) +++ branches/eam_branches/ipp-20140904/ippScripts/scripts/publish_file.pl (revision 37804) @@ -252,4 +252,5 @@ warp_id => $comp->{warp1}, diff_id => $comp->{diff_id}, + camera => $comp->{camera}, output_format => $comp->{output_format}, direction => 1, @@ -291,4 +292,6 @@ warp_id => $comp->{warp2}, diff_id => $comp->{diff_id}, + camera => $comp->{camera}, + # missing output_format? direction => 0, comment => $reg_comp->{comment}, @@ -359,20 +362,22 @@ my $name = shift; # Name of component - $data->{zp} = $comp->{zp} unless defined $data->{zp}; - $data->{zp_err} = $comp->{zp_err} unless defined $data->{zp_err}; - $data->{astrom} = $comp->{astrom} unless defined $data->{astrom}; - $data->{exp_name} = $comp->{exp_name} unless defined $data->{exp_name}; - $data->{exp_id} = $comp->{exp_id} unless defined $data->{exp_id}; - $data->{chip_id} = $comp->{chip_id} unless defined $data->{chip_id}; - $data->{cam_id} = $comp->{cam_id} unless defined $data->{cam_id}; - $data->{fake_id} = $comp->{fake_id} unless defined $data->{fake_id}; - $data->{warp_id} = $comp->{warp_id} unless defined $data->{warp_id}; - $data->{diff_id} = $comp->{diff_id} unless defined $data->{diff_id}; + $data->{zp} = $comp->{zp} unless defined $data->{zp}; + $data->{zp_err} = $comp->{zp_err} unless defined $data->{zp_err}; + $data->{astrom} = $comp->{astrom} unless defined $data->{astrom}; + $data->{exp_name} = $comp->{exp_name} unless defined $data->{exp_name}; + $data->{exp_id} = $comp->{exp_id} unless defined $data->{exp_id}; + $data->{chip_id} = $comp->{chip_id} unless defined $data->{chip_id}; + $data->{cam_id} = $comp->{cam_id} unless defined $data->{cam_id}; + $data->{fake_id} = $comp->{fake_id} unless defined $data->{fake_id}; + $data->{warp_id} = $comp->{warp_id} unless defined $data->{warp_id}; + $data->{diff_id} = $comp->{diff_id} unless defined $data->{diff_id}; + $data->{camera} = $comp->{camera} unless defined $data->{camera}; + # missing output_format? $data->{direction} = $comp->{direction} unless defined $data->{direction}; - $data->{comment} = $comp->{comment} unless defined $data->{comment}; - $data->{obsmode} = $comp->{obsmode} unless defined $data->{obsmode}; - $data->{difftype} = $comp->{difftype} unless defined $data->{difftype}; - $data->{sky} = $comp->{sky} unless defined $data->{sky}; - $data->{shutoutc} = $comp->{shutoutc} unless defined $data->{shutoutc}; + $data->{comment} = $comp->{comment} unless defined $data->{comment}; + $data->{obsmode} = $comp->{obsmode} unless defined $data->{obsmode}; + $data->{difftype} = $comp->{difftype} unless defined $data->{difftype}; + $data->{sky} = $comp->{sky} unless defined $data->{sky}; + $data->{shutoutc} = $comp->{shutoutc} unless defined $data->{shutoutc}; &my_die("zp value for $name doesn't match", $pub_id, $PS_EXIT_SYS_ERROR) if defined $data->{zp} and $comp->{zp} != $data->{zp}; @@ -386,4 +391,5 @@ &my_die("warp_id value for $name doesn't match", $pub_id, $PS_EXIT_SYS_ERROR) if defined $data->{warp_id} and $comp->{warp_id} != $data->{warp_id}; &my_die("diff_id value for $name doesn't match", $pub_id, $PS_EXIT_SYS_ERROR) if defined $data->{diff_id} and $comp->{diff_id} != $data->{diff_id}; + &my_die("camera value for $name doesn't match", $pub_id, $PS_EXIT_SYS_ERROR) if defined $data->{camera} and $comp->{camera} ne $data->{camera}; &my_die("direction value for $name doesn't match", $pub_id, $PS_EXIT_SYS_ERROR) if defined $data->{direction} and $comp->{direction} != $data->{direction}; &my_die("comment value for $name doesn't match", $pub_id, $PS_EXIT_SYS_ERROR) if defined $data->{comment} and $comp->{comment} ne $data->{comment}; @@ -414,4 +420,5 @@ $command .= " -warp_id " . $data->{warp_id} if defined $data->{warp_id}; $command .= " -diff_id " . $data->{diff_id} if defined $data->{diff_id}; + $command .= " -camera " . $data->{camera} if defined $data->{camera}; $command .= " -inverse" if defined $data->{direction} and $data->{direction} == 0; $command .= " -zp " . $data->{zp} if defined $data->{zp}; @@ -420,5 +427,5 @@ $command .= " -version " . $data->{output_format} if defined $data->{output_format}; $command .= " -comment \"" . $data->{comment} . "\"" if defined $data->{comment}; - $command .= " -obsmode " . $data->{obsmode} if defined $data->{obsmode}; + $command .= " -obsmode \"" . $data->{obsmode} . "\"" if defined $data->{obsmode}; $command .= " -difftype " . $data->{difftype} if defined $data->{difftype}; $command .= " -sky " . $data->{sky} if defined $data->{sky}; Index: branches/eam_branches/ipp-20140904/ippScripts/scripts/rawcheck.pl =================================================================== --- branches/eam_branches/ipp-20140904/ippScripts/scripts/rawcheck.pl (revision 37803) +++ branches/eam_branches/ipp-20140904/ippScripts/scripts/rawcheck.pl (revision 37804) @@ -62,8 +62,11 @@ my %backup_hosts = ('ippb00' => 1, 'ippb01' => 1, 'ippb02' => 1, 'ippb03' => 1, - 'ippb04' => 1, 'ippb05' => 1 + 'ippb04' => 1, 'ippb05' => 1, + 'ippb06' => 1 ); -my %backup_destinations = ('ippb04' => 1, - 'ippb05' => 1); +my %backup_destinations = (# 'ippb04' => 1, # full + # 'ippb05' => 1, # full + 'ippb06' => 1 + ); my $backup_Nvols = 3; Index: branches/eam_branches/ipp-20140904/ippScripts/scripts/stack_bkg_mk_mdc.pl =================================================================== --- branches/eam_branches/ipp-20140904/ippScripts/scripts/stack_bkg_mk_mdc.pl (revision 37803) +++ branches/eam_branches/ipp-20140904/ippScripts/scripts/stack_bkg_mk_mdc.pl (revision 37804) @@ -13,4 +13,8 @@ my ($dbname,$stack_id,$label,$camera,$filter); + +my $ipprc = PS::IPP::Config->new(); # IPP Configuration +my $siteConfig = $ipprc->{_siteConfig}; + #$dbname = 'gpc1'; #$camera = 'GPC1'; @@ -29,7 +33,10 @@ use constant DB_SOCKET => '/var/run/mysqld/mysqld.sock'; -my $dbserver = 'ippdb01'; -my $dbuser = 'ippuser'; -my $dbpass = 'ippuser'; +my $dbserver = metadataLookupStr($siteConfig, 'DBSERVER'); +my $dbuser = metadataLookupStr($ipprc->{_siteConfig}, "RO_DBUSER"); +my $dbpass = metadataLookupStr($ipprc->{_siteConfig}, "RO_DBPASSWORD"); +die "database configuration not set up" unless defined($dbserver); +die "database configuration not set up" unless defined($dbuser); +die "database configuration not set up" unless defined($dbpass); my $db = DBI->connect("DBI:mysql:database=${dbname};host=${dbserver};" . "mysql_socket=" . DB_SOCKET(), @@ -38,5 +45,5 @@ ) or die "Unable to connect to database $DBI::errstr\n"; -my $ipprc = PS::IPP::Config->new( $camera ); +$ipprc = PS::IPP::Config->new( $camera ); my $mdcParser = PS::IPP::Metadata::Config->new; Index: branches/eam_branches/ipp-20140904/ippTasks/addstar.pro =================================================================== --- branches/eam_branches/ipp-20140904/ippTasks/addstar.pro (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTasks/addstar.pro (revision 37804) @@ -52,4 +52,10 @@ active false end + task addstar.revert.diff + active false + end + task addstar.revert.fullforce + active false + end end @@ -65,4 +71,10 @@ end task addstar.revert.skycal + active true + end + task addstar.revert.diff + active true + end + task addstar.revert.fullforce active true end @@ -141,4 +153,6 @@ $addstar_revert_DB_SSM = 0 $addstar_revert_DB_SC = 0 +$addstar_revert_DB_FF = 0 +$addstar_revert_DB_DF = 0 #$LOADEXEC_ADD = 300 @@ -275,4 +289,10 @@ sprintf outroot "%s.add.%s" $STAGEROOT $ADD_ID end + if ("$STAGE" == "diff") + sprintf outroot "%s.add.%s" $STAGEROOT $ADD_ID + end + if ("$STAGE" == "fullforce") + sprintf outroot "%s.add.%s" $STAGEROOT $ADD_ID + end @@ -293,4 +313,10 @@ if ("$STAGE" == "skycal") $run = $run --stage_id $STAGE_ID + end + if ("$STAGE" == "diff") + $run = $run --stage_id $STAGE_ID --stage_extra1 $STAGE_EXTRA1 + end + if ("$STAGE" == "fullforce") + $run = $run --stage_id $STAGE_ID --stage_extra1 $STAGE_EXTRA1 end @@ -526,2 +552,98 @@ end end + + +task addstar.revert.diff + host local + + periods -poll 5.0 + periods -exec 60.0 + periods -timeout 1200 + npending 1 + active false + + stdout NULL + stderr $LOGDIR/revert.log + + task.exec + if ($LABEL:n == 0) break + $run = addtool -revertprocessedexp -stage diff + if ($DB:n == 0) + option DEFAULT + else + # save the DB name for the exit tasks + option $DB:$addstar_revert_DB_DF + $run = $run -dbname $DB:$addstar_revert_DB_DF + $addstar_revert_DB_DF ++ + if ($addstar_revert_DB_DF >= $DB:n) set addstar_revert_DB_DF = 0 + end + add_poll_labels run + command $run + end + + # success + task.exit 0 + end + + # locked list + task.exit default + showcommand failure + end + + task.exit crash + showcommand crash + end + + # operation times out? + task.exit timeout + showcommand timeout + end +end + + +task addstar.revert.fullforce + host local + + periods -poll 5.0 + periods -exec 60.0 + periods -timeout 1200 + npending 1 + active false + + stdout NULL + stderr $LOGDIR/revert.log + + task.exec + if ($LABEL:n == 0) break + $run = addtool -revertprocessedexp -stage fullforce + if ($DB:n == 0) + option DEFAULT + else + # save the DB name for the exit tasks + option $DB:$addstar_revert_DB_FF + $run = $run -dbname $DB:$addstar_revert_DB_FF + $addstar_revert_DB_FF ++ + if ($addstar_revert_DB_FF >= $DB:n) set addstar_revert_DB_FF = 0 + end + add_poll_labels run + command $run + end + + # success + task.exit 0 + end + + # locked list + task.exit default + showcommand failure + end + + task.exit crash + showcommand crash + end + + # operation times out? + task.exit timeout + showcommand timeout + end +end Index: branches/eam_branches/ipp-20140904/ippTasks/register.pro =================================================================== --- branches/eam_branches/ipp-20140904/ippTasks/register.pro (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTasks/register.pro (revision 37804) @@ -36,4 +36,7 @@ active false end + task register.exp.revert + active false + end task register.exp.load active true @@ -60,4 +63,7 @@ active false end + task register.exp.revert + active false + end task register.exp.load active false @@ -78,8 +84,14 @@ active true end + task register.exp.revert + active true + end end macro register.revert.off task register.imfile.revert + active false + end + task register.exp.revert active false end @@ -202,4 +214,5 @@ $regPendingBurntoolImfile_DB = 0 $regRevertImfile_DB = 0 +$regRevertExp_DB = 0 $regPendingExp_DB = 0 $reg_datePAGE = 0 @@ -371,4 +384,50 @@ $regRevertImfile_DB ++ if ($regRevertImfile_DB >= $DB:n) set regRevertImfile_DB = 0 + end + add_poll_args run + command $run + end + + # success + task.exit 0 + end + + # locked list + task.exit default + showcommand failure + end + + task.exit crash + showcommand crash + end + + # operation times out + task.exit timeout + showcommand timeout + end +end + + + +task register.exp.revert + host local + periods -poll 60.0 + periods -exec 1800.0 + periods -timeout 120.0 + npending 1 + + stdout NULL + stderr $LOGDIR/register.exp.revert.log + + task.exec +# Only revert fault 2, because other faults likely mean something horrible is happening. + $run = regtool -revertprocessedexp -fault 2 + if ($DB:n == 0) + option DEFAULT + else + option $DB:$regRevertExp_DB + $run = $run -dbname $DB:$regRevertExp_DB + $regRevertExp_DB ++ + if ($regRevertExp_DB >= $DB:n) set regRevertExp_DB = 0 end add_poll_args run Index: branches/eam_branches/ipp-20140904/ippTasks/stack.pro =================================================================== --- branches/eam_branches/ipp-20140904/ippTasks/stack.pro (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTasks/stack.pro (revision 37804) @@ -53,10 +53,4 @@ active false end - task stack.summary.load - active true - end - task stack.summary.run - active true - end end @@ -72,12 +66,7 @@ active false end - task stack.summary.load - active false - end - task stack.summary.run - active false - end -end - +end + +# Same for summary stages. macro stack.summary.on task stack.summary.load @@ -313,5 +302,6 @@ periods -exec $LOADEXEC periods -timeout 1200 - trange 07:00:00 08:00:00 -nmax 1 + active false +# trange 07:00:00 08:00:00 -nmax 1 npending 1 @@ -368,4 +358,5 @@ periods -exec $RUNEXEC periods -timeout 60 + active false task.exec Index: branches/eam_branches/ipp-20140904/ippTasks/survey.pro =================================================================== --- branches/eam_branches/ipp-20140904/ippTasks/survey.pro (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTasks/survey.pro (revision 37804) @@ -178,13 +178,14 @@ # user functions to manipulate warp/stack diff labels macro survey.add.WSdiff - if ($0 != 6) - echo "USAGE: survey.add.WSdiff (warp label) (stack label) (dist_group) (workdir base) (target label)" + if ($0 != 7) + echo "USAGE: survey.add.WSdiff (key) (warp label) (stack label) (dist_group) (workdir base) (target label)" break end book newpage SURVEY_DIFF_WARPSTACK $1 - book setword SURVEY_DIFF_WARPSTACK $1 STACK_LABEL $2 - book setword SURVEY_DIFF_WARPSTACK $1 DIST_GROUP $3 - book setword SURVEY_DIFF_WARPSTACK $1 WORKDIR $4 - book setword SURVEY_DIFF_WARPSTACK $1 TARGET_LABEL $5 + book setword SURVEY_DIFF_WARPSTACK $1 WARP_LABEL $2 + book setword SURVEY_DIFF_WARPSTACK $1 STACK_LABEL $3 + book setword SURVEY_DIFF_WARPSTACK $1 DIST_GROUP $4 + book setword SURVEY_DIFF_WARPSTACK $1 WORKDIR $5 + book setword SURVEY_DIFF_WARPSTACK $1 TARGET_LABEL $6 book setword SURVEY_DIFF_WARPSTACK $1 STATE PENDING end @@ -697,4 +698,7 @@ periods -exec $SURVEY_EXEC periods -timeout $SURVEY_TIMEOUT + ## time in UT -- want to start ~0900 HST and run until next night starts. 24h boundary so two trange needed? + trange 19:00:00 23:59:00 + trange 00:00:00 05:00:00 npending 1 @@ -710,13 +714,13 @@ endif - book getpage SURVEY_DIFF_WARPSTACK 0 -var label -key STATE NEW - if ("$label" == "NULL") + book getpage SURVEY_DIFF_WARPSTACK 0 -var wsKey -key STATE NEW + if ("$wsKey" == "NULL") # All labels have been done --- reset # echo "Resetting labels" for i 0 $N - book getpage SURVEY_DIFF_WARPSTACK $i -var label - book setword SURVEY_DIFF_WARPSTACK $label STATE NEW + book getpage SURVEY_DIFF_WARPSTACK $i -var wsKey + book setword SURVEY_DIFF_WARPSTACK $wsKey STATE NEW end - book getpage SURVEY_DIFF_WARPSTACK 0 -var label -key STATE NEW + book getpage SURVEY_DIFF_WARPSTACK 0 -var wsKey -key STATE NEW # Select different database @@ -725,9 +729,10 @@ end - book setword SURVEY_DIFF_WARPSTACK $label STATE DONE - book getword SURVEY_DIFF_WARPSTACK $label STACK_LABEL -var stack_label - book getword SURVEY_DIFF_WARPSTACK $label WORKDIR -var workdir - book getword SURVEY_DIFF_WARPSTACK $label DIST_GROUP -var dist_group - book getword SURVEY_DIFF_WARPSTACK $label TARGET_LABEL -var target_label + book setword SURVEY_DIFF_WARPSTACK $wsKey STATE DONE + book getword SURVEY_DIFF_WARPSTACK $wsKey WARP_LABEL -var warp_label + book getword SURVEY_DIFF_WARPSTACK $wsKey STACK_LABEL -var stack_label + book getword SURVEY_DIFF_WARPSTACK $wsKey WORKDIR -var workdir + book getword SURVEY_DIFF_WARPSTACK $wsKey DIST_GROUP -var dist_group + book getword SURVEY_DIFF_WARPSTACK $wsKey TARGET_LABEL -var target_label $year = `date +%Y` @@ -737,5 +742,5 @@ # Do we need anything else here? $run = difftool -definewarpstack -good_frac 0.2 - $run = $run -warp_label $label + $run = $run -warp_label $warp_label $run = $run -stack_label $stack_label $run = $run -set_dist_group $dist_group @@ -1012,5 +1017,10 @@ $run = $run -uncensored end - + if ("$stage" == "fullforce") + $run = $run -uncensored + end + if ("$stage" == "diff") + $run = $run -uncensored) + end # we need to handle isp/gpc1: gpc1 wants destreaked, isp wants uncensored. Perhaps the best way is by default do destreaked (the most paranoid), and if db = isp (is this the best way?) then do uncensored Index: branches/eam_branches/ipp-20140904/ippTools/share/Makefile.am =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/Makefile.am (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTools/share/Makefile.am (revision 37804) @@ -8,4 +8,8 @@ addtool_find_cam_id_dvo.sql \ addtool_find_cam_id.sql \ + addtool_find_diff_id.sql \ + addtool_find_diff_id_dvo.sql \ + addtool_find_ff_id.sql \ + addtool_find_ff_id_dvo.sql \ addtool_find_skycal_id_dvo.sql \ addtool_find_skycal_id.sql \ @@ -16,4 +20,6 @@ addtool_find_pendingexp.sql \ addtool_find_pendingexp_cam.sql \ + addtool_find_pendingexp_diff.sql \ + addtool_find_pendingexp_ff.sql \ addtool_find_pendingexp_skycal.sql \ addtool_find_pendingexp_stack.sql \ @@ -21,4 +27,6 @@ addtool_find_pendingmergeprocess.sql \ addtool_find_processedexp_cam.sql \ + addtool_find_processedexp_diff.sql \ + addtool_find_processedexp_ff.sql \ addtool_find_processedexp_skycal.sql \ addtool_find_processedexp_stack.sql \ @@ -29,4 +37,6 @@ addtool_pendingcleanuprun.sql \ addtool_queue_cam_id.sql \ + addtool_queue_diff_id.sql \ + addtool_queue_ff_id.sql \ addtool_queue_stack_id.sql \ addtool_queue_sky_id_multi.sql \ @@ -35,4 +45,6 @@ addtool_revertminidvodbprocessed.sql \ addtool_revertprocessedexp_cam.sql \ + addtool_revertprocessedexp_diff.sql \ + addtool_revertprocessedexp_ff.sql \ addtool_revertprocessedexp_skycal.sql \ addtool_revertprocessedexp_stack.sql \ @@ -471,4 +483,5 @@ laptool_revertgroup.sql \ laptool_filtersforgroup.sql \ + laptool_revertsass.sql \ vptool_find_rawexp.sql \ vptool_pendingimfile.sql \ Index: branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_diff_id.sql =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_diff_id.sql (revision 37804) +++ branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_diff_id.sql (revision 37804) @@ -0,0 +1,22 @@ +SELECT + diffRun.*, diffInputSkyfile.diff_skyfile_id +FROM + diffInputSkyfile +JOIN diffSkyfile using (diff_id, skycell_id) +JOIN diffRun using (diff_id) + +LEFT JOIN (SELECT diff_id AS added_diff_id, + addRun.dvodb AS previous_dvodb + FROM addRun +JOIN diffRun on (diff_id = stage_id and addRun.stage = 'diff') + ) as foo + ON diff_id = added_diff_id + -- AND stage = 'diff' + -- hook for qualifying the join on the previous_dvodb + AND %s +WHERE + diffRun.state = 'full' + AND diffSkyfile.fault = 0 + AND diffSkyfile.quality = 0 + AND added_diff_id IS NULL + -- addtool adds checks on exposure being added to the dvodb previously Index: branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_diff_id_dvo.sql =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_diff_id_dvo.sql (revision 37804) +++ branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_diff_id_dvo.sql (revision 37804) @@ -0,0 +1,12 @@ +SELECT diffRun.*, diffInputSkyfile.diff_skyfile_id +FROM + diffInputSkyfile +JOIN diffSkyfile using (diff_id, skycell_id) +JOIN diffRun using (diff_id) + +WHERE diffRun.state = 'full' and diffSkyfile.quality = 0 + AND diffRun.diff_id NOT IN (SELECT diffRun.diff_id + FROM addRun + JOIN diffRun on (diffRun.diff_id = addRun.stage_id and addRun.stage = 'diff') JOIN diffSkyfile on (diffRun.diff_id = diffSkyfile.diff_id) + WHERE addRun.stage = 'diff' AND %s + ) Index: branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_ff_id.sql =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_ff_id.sql (revision 37804) +++ branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_ff_id.sql (revision 37804) @@ -0,0 +1,21 @@ +SELECT + fullForceRun.*, fullForceResult.warp_id +FROM + fullForceRun +JOIN fullForceResult using (ff_id) + +LEFT JOIN (SELECT ff_id AS added_ff_id, + addRun.dvodb AS previous_dvodb + FROM addRun +JOIN fullForceRun on (ff_id = stage_id and addRun.stage = 'fullforce') + ) as foo + ON ff_id = added_ff_id + --AND stage = 'skycal' + -- hook for qualifying the join on the previous_dvodb + AND %s +WHERE + fullForceRun.state = 'full' + AND fullForceResult.quality = 0 + AND fullForceResult.fault = 0 + AND added_ff_id IS NULL + -- addtool adds checks on exposure being added to the dvodb previously Index: branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_ff_id_dvo.sql =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_ff_id_dvo.sql (revision 37804) +++ branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_ff_id_dvo.sql (revision 37804) @@ -0,0 +1,8 @@ +SELECT fullForceRun.*, fullForceResult.warp_id FROM fullForceRun JOIN fullForceResult using (ff_id) + +WHERE fullForceRun.state = 'full' and fullForceResult.quality = 0 + AND fullForceRun.ff_id NOT IN (SELECT fullForceRun.ff_id + FROM addRun + JOIN fullForceRun on (fullForceRun.ff_id = addRun.stage_id and addRun.stage = 'fullforce') JOIN fullForceResult on (fullForceRun.ff_id = fullForceResult.ff_id) + WHERE addRun.stage = 'fullforce' AND %s + ) Index: branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_pendingexp_diff.sql =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_pendingexp_diff.sql (revision 37804) +++ branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_pendingexp_diff.sql (revision 37804) @@ -0,0 +1,21 @@ + SELECT + addRun.*, + diffSkyfile.path_base as stageroot, + 'GPC1' as camera, + 'PS1' as telescope +FROM addRun +JOIN diffInputSkyfile + ON (diff_id = stage_id and diff_skyfile_id = stage_extra1 and stage = "diff") +JOIN diffSkyfile + USING (diff_id, skycell_id) +JOIN diffRun + using(diff_id) +LEFT JOIN addProcessedExp using (add_id) +LEFT JOIN addMask + ON addRun.label = addMask.label +WHERE + diffRun.state = 'full' + AND ((addRun.state = 'new' AND addProcessedExp.add_id IS NULL) OR addRun.state = 'update') + AND addRun.dvodb IS NOT NULL + AND addRun.workdir IS NOT NULL + AND addMask.label IS NULL Index: branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_pendingexp_ff.sql =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_pendingexp_ff.sql (revision 37804) +++ branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_pendingexp_ff.sql (revision 37804) @@ -0,0 +1,19 @@ + SELECT + addRun.*, + fullForceResult.path_base as stageroot, + 'GPC1' as camera, + 'PS1' as telescope +FROM addRun +JOIN fullForceResult + ON (ff_id = stage_id and warp_id = stage_extra1 and stage = "fullforce") +JOIN fullForceRun + using(ff_id) +LEFT JOIN addProcessedExp using (add_id) +LEFT JOIN addMask + ON addRun.label = addMask.label +WHERE + fullForceRun.state = 'full' + AND ((addRun.state = 'new' AND addProcessedExp.add_id IS NULL) OR addRun.state = 'update') + AND addRun.dvodb IS NOT NULL + AND addRun.workdir IS NOT NULL + AND addMask.label IS NULL Index: branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_processedexp_diff.sql =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_processedexp_diff.sql (revision 37804) +++ branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_processedexp_diff.sql (revision 37804) @@ -0,0 +1,6 @@ +SELECT + addProcessedExp.*, + addRun.workdir +FROM addProcessedExp +JOIN addRun + USING(add_id) Index: branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_processedexp_ff.sql =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_processedexp_ff.sql (revision 37804) +++ branches/eam_branches/ipp-20140904/ippTools/share/addtool_find_processedexp_ff.sql (revision 37804) @@ -0,0 +1,6 @@ +SELECT + addProcessedExp.*, + addRun.workdir +FROM addProcessedExp +JOIN addRun + USING(add_id) Index: branches/eam_branches/ipp-20140904/ippTools/share/addtool_queue_diff_id.sql =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/addtool_queue_diff_id.sql (revision 37804) +++ branches/eam_branches/ipp-20140904/ippTools/share/addtool_queue_diff_id.sql (revision 37804) @@ -0,0 +1,23 @@ +INSERT INTO addRun + SELECT + 0, -- add_id + 'diff', -- stage + diff_id, -- stage_id + %d, -- stage_extra1 + '%s', -- state + '%s', -- workdir + '%s', -- workdir_state + '%s', -- reduction + '%s', -- label + '%s', -- data_group + '%s', -- dvodb + '%s', -- note + %d, -- image_only + %d, -- minidvodb + '%s', -- minidvodb_group + '%s' -- minidvodb_name + FROM diffRun + WHERE + diffRun.state = 'full' + AND diffRun.diff_id = %lld + Index: branches/eam_branches/ipp-20140904/ippTools/share/addtool_queue_ff_id.sql =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/addtool_queue_ff_id.sql (revision 37804) +++ branches/eam_branches/ipp-20140904/ippTools/share/addtool_queue_ff_id.sql (revision 37804) @@ -0,0 +1,22 @@ +INSERT INTO addRun + SELECT + 0, -- add_id + 'fullforce', -- stage + ff_id, -- stage_id + %d, -- stage_extra1 + '%s', -- state + '%s', -- workdir + '%s', -- workdir_state + '%s', -- reduction + '%s', -- label + '%s', -- data_group + '%s', -- dvodb + '%s', -- note + %d, -- image_only + %d, -- minidvodb + '%s', -- minidvodb_group + '%s' -- minidvodb_name + FROM fullForceRun + WHERE + fullForceRun.state = 'full' + AND fullForceRun.ff_id = %lld Index: branches/eam_branches/ipp-20140904/ippTools/share/addtool_revertprocessedexp_diff.sql =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/addtool_revertprocessedexp_diff.sql (revision 37804) +++ branches/eam_branches/ipp-20140904/ippTools/share/addtool_revertprocessedexp_diff.sql (revision 37804) @@ -0,0 +1,7 @@ +DELETE FROM addProcessedExp +USING addProcessedExp, addRun +WHERE + addRun.add_id = addProcessedExp.add_id + AND addProcessedExp.fault != 0 + AND addRun.state = 'new' + AND addRun.stage = 'diff' Index: branches/eam_branches/ipp-20140904/ippTools/share/addtool_revertprocessedexp_ff.sql =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/addtool_revertprocessedexp_ff.sql (revision 37804) +++ branches/eam_branches/ipp-20140904/ippTools/share/addtool_revertprocessedexp_ff.sql (revision 37804) @@ -0,0 +1,7 @@ +DELETE FROM addProcessedExp +USING addProcessedExp, addRun +WHERE + addRun.add_id = addProcessedExp.add_id + AND addProcessedExp.fault != 0 + AND addRun.state = 'new' + AND addRun.stage = 'fullforce' Index: branches/eam_branches/ipp-20140904/ippTools/share/fftool_result.sql =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/fftool_result.sql (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTools/share/fftool_result.sql (revision 37804) @@ -11,4 +11,6 @@ stackRun.skycell_id, stackRun.filter, + rawExp.exp_id, + rawExp.exp_name, rawExp.camera FROM fullForceRun Index: branches/eam_branches/ipp-20140904/ippTools/share/laptool_exposures.sql =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/laptool_exposures.sql (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTools/share/laptool_exposures.sql (revision 37804) @@ -1,7 +1,10 @@ -select DISTINCT V3.*, - diffRun.diff_id,diffRun.state as diffRun_state, - coalesce(CONVERT(sum(diffSkyfile.quality != 0),SIGNED),0) AS diff_bad_quality, - coalesce(CONVERT(count(diffSkyfile.diff_id),SIGNED),0) AS diff_component_count - FROM +-- Join the diff information +-- select DISTINCT V3.*, +-- diffRun.diff_id,diffRun.state as diffRun_state, +-- coalesce(CONVERT(sum(diffSkyfile.quality != 0),SIGNED),0) AS diff_bad_quality, +-- coalesce(CONVERT(count(diffSkyfile.diff_id),SIGNED),0) AS diff_component_count +-- FROM +select V3.*,-1 as diff_id ,"full" as diffRun_state,0 as diff_bad_quality,0 AS diff_component_count FROM +-- Join the warp information ( select V2.*, warpRun.warp_id,warpRun.state as warpRun_state, @@ -10,4 +13,5 @@ warpRun.magicked FROM +-- Join the camera and fake information ( select V1.*, camRun.cam_id,camRun.state as camRun_state, @@ -15,4 +19,5 @@ coalesce(CONVERT(count(camProcessedExp.cam_id),SIGNED),0) AS cam_component_count, fakeRun.fake_id,fakeRun.state as fakeRun_state FROM +-- Get all the lap, exposure, chip stage information ( SELECT DISTINCT lap_id,lapRun.tess_id,projection_cell,filter,lapRun.state as lapRun_state, lapRun.registered, lapRun.fault, lapRun.label, lapRun.dist_group, @@ -26,16 +31,19 @@ GROUP BY lap_id,exp_id ) AS V1 +-- End lap/exposure/chip LEFT JOIN camRun USING(chip_id) LEFT JOIN camProcessedExp USING(cam_id) LEFT JOIN fakeRun USING(cam_id) GROUP BY lap_id,exp_id ) AS V2 +-- End camera/fake LEFT JOIN warpRun USING(fake_id) LEFT JOIN warpSkyfile USING(warp_id) GROUP BY lap_id,exp_id ) AS V3 -LEFT JOIN - (SELECT DISTINCT diff_id,warp1,warp2 FROM diffInputSkyfile) AS DI ON - (DI.warp1 = warp_id OR DI.warp2 = warp_id) -LEFT JOIN diffRun USING(diff_id) -LEFT JOIN diffSkyfile USING(diff_id) +-- End warp +-- LEFT JOIN +-- (SELECT DISTINCT diff_id,warp1,warp2 FROM diffInputSkyfile) AS DI ON +-- (DI.warp1 = warp_id OR DI.warp2 = warp_id) +-- LEFT JOIN diffRun USING(diff_id) +-- LEFT JOIN diffSkyfile USING(diff_id) GROUP BY lap_id,exp_id Index: branches/eam_branches/ipp-20140904/ippTools/share/laptool_revertsass.sql =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/laptool_revertsass.sql (revision 37804) +++ branches/eam_branches/ipp-20140904/ippTools/share/laptool_revertsass.sql (revision 37804) @@ -0,0 +1,1 @@ +UPDATE lapRun SET final_sass_id = NULL,quick_sass_id = NULL Index: branches/eam_branches/ipp-20140904/ippTools/share/pstamptool_pendingdependent.sql =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/share/pstamptool_pendingdependent.sql (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTools/share/pstamptool_pendingdependent.sql (revision 37804) @@ -2,4 +2,5 @@ pstampDependent.*, pstampRequest.label, + MIN(req_id) as first_req_id, concat( pstampDependent.outdir, '/checkdep.', dep_id, '.log') as logfile, IFNULL(Label.priority, 10000) AS priority Index: branches/eam_branches/ipp-20140904/ippTools/src/addtool.c =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/src/addtool.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTools/src/addtool.c (revision 37804) @@ -122,4 +122,7 @@ PXOPT_COPY_S64(config->args, where, "-sky_id", "staticskyRun.sky_id", "=="); PXOPT_COPY_S64(config->args, where, "-skycal_id", "skycalRun.skycal_id", "=="); + PXOPT_COPY_S64(config->args, where, "-diff_id", "diffRun.diff_id", "=="); + PXOPT_COPY_S64(config->args, where, "-ff_id", "fullForceRun.ff_id", "=="); + PXOPT_LOOKUP_STR(stage, config->args, "-stage", false, false); if (strcmp(stage, "cam")== 0) { @@ -146,5 +149,17 @@ PXOPT_COPY_STR(config->args, where, "-reduction", "skycalRun.reduction", "=="); } - + if (strcmp(stage, "diff")==0) { + pxAddLabelSearchArgs (config, where, "-label", "diffRun.label", "=="); //define using skycalRun label + pxAddLabelSearchArgs (config, where, "-data_group","diffRun.data_group", "=="); + PXOPT_COPY_STR(config->args, where, "-reduction", "diffRun.reduction", "=="); + } + if (strcmp(stage, "fullforce")==0) { + pxAddLabelSearchArgs (config, where, "-label", "fullForceRun.label", "=="); //define using skycalRun label + pxAddLabelSearchArgs (config, where, "-data_group","fullForceRun.data_group", "=="); + PXOPT_COPY_STR(config->args, where, "-reduction", "fullForceRun.reduction", "=="); + } + + + if (!psListLength(where->list)) { psFree(where); @@ -258,5 +273,37 @@ } - + if (strcmp(stage,"diff") == 0) { + if (dvodb ) { + psTrace("addtool.c", PS_LOG_INFO, "dvodb argument found (%s) using addtool_find_diff_id_dvo.sql\n%s\n", dvodb,stage); + // find the skycal_id of all the exposures that we want to queue up. + bare_query = pxDataGet("addtool_find_diff_id_dvo.sql"); + // user supplied dvodb + psStringAppend(&dvodb_string, "addRun.dvodb = '%s'", dvodb); + } else { + psTrace("addtool.c", PS_LOG_INFO, "dvodb argument not found using addtool_find_diff_id.sql\n%s\n",stage); + // find the skycal_id of all the exposures that we want to queue up. + bare_query = pxDataGet("addtool_find_diff_id.sql"); + // inherit dvodb from skycalRun, avoid matching NULL + psStringAppend(&dvodb_string, "(diffRun.dvodb IS NOT NULL AND previous_dvodb = diffRun.dvodb)"); + // this is silly, there is no dvodb in skycalRun...? + } + } + + if (strcmp(stage,"fullforce") == 0) { + if (dvodb ) { + psTrace("addtool.c", PS_LOG_INFO, "dvodb argument found (%s) using addtool_find_ff_id_dvo.sql\n%s\n", dvodb,stage); + // find the skycal_id of all the exposures that we want to queue up. + bare_query = pxDataGet("addtool_find_ff_id_dvo.sql"); + // user supplied dvodb + psStringAppend(&dvodb_string, "addRun.dvodb = '%s'", dvodb); + } else { + psTrace("addtool.c", PS_LOG_INFO, "dvodb argument not found using addtool_ff_id.sql\n%s\n",stage); + // find the skycal_id of all the exposures that we want to queue up. + bare_query = pxDataGet("addtool_find_ff_id.sql"); + // inherit dvodb from skycalRun, avoid matching NULL + psStringAppend(&dvodb_string, "(fullForceRun.dvodb IS NOT NULL AND previous_dvodb = fullForceRun.dvodb)"); + // this is silly, there is no dvodb in skycalRun...? + } + } @@ -324,4 +371,15 @@ psStringAppend(&query, " GROUP BY skycal_id, sky_id, stack_id "); //this needs checking, but I think it shoul be fine? it groups by lots of stuff (including stack - we only want one of each stack in there } + //needs to be checked HAF xxx + if (strcmp(stage,"diff") == 0) { + psStringAppend(&query, " GROUP BY diff_id, diff_skyfile_id "); //this needs checking, but I think it shoul be fine? it groups by lots of stuff (including stack - we only want one of each stack in there + } + //needs to be checked HAF xxx + if (strcmp(stage,"fullforce") == 0) { + psStringAppend(&query, " GROUP BY ff_id, warp_id "); //this needs checking, but I think it shoul be fine? it groups by lots of stuff (including stack - we only want one of each stack in there + } + + + psTrace("addtool.c", PS_LOG_INFO,"query: \n\n%s\n\n",query); @@ -476,4 +534,58 @@ } + if (strcmp(stage,"diff") == 0) { + + for (long i = 0; i < psArrayLength(output); i++) { + psMetadata *md = output->data[i]; + + diffRunRow *row = diffRunObjectFromMetadata(md); + + if (!row) { + psError(PS_ERR_UNKNOWN, false, "failed to convert metadata into diffRun"); + psFree(output); + return false; + } + + if (!dvodb) { // there's no skycalRun.dvodb + psError(PS_ERR_UNKNOWN, false, "cannot queue addstar run without a defined dvodb: label: %s, diff_id %" PRId64, row->label, row->diff_id); + psFree(output); + return false; + } + if (!workdir && !row->workdir) { + psError(PS_ERR_UNKNOWN, false, "cannot queue addstar run without a defined workdir: label: %s, diff_id %" PRId64, row->label, row->diff_id); + psFree(output); + return false; + } + + psFree(row); + } + } + + if (strcmp(stage,"fullforce") == 0) { + for (long i = 0; i < psArrayLength(output); i++) { + psMetadata *md = output->data[i]; + + fullForceRunRow *row = fullForceRunObjectFromMetadata(md); + + if (!row) { + psError(PS_ERR_UNKNOWN, false, "failed to convert metadata into fullforceRun"); + psFree(output); + return false; + } + + if (!dvodb) { // there's no skycalRun.dvodb + psError(PS_ERR_UNKNOWN, false, "cannot queue addstar run without a defined dvodb: label: %s, ff_id %" PRId64, row->label, row->ff_id); + psFree(output); + return false; + } + if (!workdir && !row->workdir) { + psError(PS_ERR_UNKNOWN, false, "cannot queue addstar run without a defined workdir: label: %s, ff_id %" PRId64, row->label, row->ff_id); + psFree(output); + return false; + } + + psFree(row); + } + } @@ -672,4 +784,106 @@ } + if (strcmp(stage,"diff") == 0) { + + for (long i = 0; i < psArrayLength(output); i++) { + psMetadata *md = output->data[i]; + psS64 stage_id =0; + // psS64 stage_extra1 = 0; + diffRunRow *row = diffRunObjectFromMetadata(md); + stage_id = row->diff_id; + bool status = false; + psS32 stage_extra1 = psMetadataLookupS32(&status, md, "diff_skyfile_id"); + if (!status) { + psError(PS_ERR_UNKNOWN, true, "failed to lookup value for item diff_skyfile_id"); + return false; + } + + + if (!row) { + psError(PS_ERR_UNKNOWN, false, "failed to convert metadata into skycalRun"); + psFree(output); + return false; + } + + // queue the exp + if (!pxaddQueueByCamID(config, + stage, + stage_id, + stage_extra1, + workdir ? workdir : row->workdir, + reduction ? reduction : row->reduction, + label ? label : row->label, + data_group ? data_group : (row->data_group ? row->data_group : (label ? label : row->label)), + dvodb ? dvodb : NULL, + note ? note : NULL, + image_only, + minidvodb, + minidvodb_group, + minidvodb_name + )) { + if (!psDBRollback(config->dbh)) { + psError(PS_ERR_UNKNOWN, false, "database error sfg"); + } + psError(PS_ERR_UNKNOWN, false, + "failed to trying to queue stage %s %" PRId64,stage, stage_id); + psFree(row); + psFree(output); + return false; + } + psFree(row); + } + } + + if (strcmp(stage,"fullforce") == 0) { + for (long i = 0; i < psArrayLength(output); i++) { + psMetadata *md = output->data[i]; + psS64 stage_id =0; + + fullForceRunRow *row = fullForceRunObjectFromMetadata(md); + stage_id = row->ff_id; + bool status = false; + psS32 stage_extra1 = psMetadataLookupS32(&status, md, "warp_id"); + if (!status) { + psError(PS_ERR_UNKNOWN, true, "failed to lookup value for item warp_id"); + return false; + } + + if (!row) { + psError(PS_ERR_UNKNOWN, false, "failed to convert metadata into skycalRun"); + psFree(output); + return false; + } + + // queue the exp + if (!pxaddQueueByCamID(config, + stage, + stage_id, + stage_extra1, + workdir ? workdir : row->workdir, + reduction ? reduction : row->reduction, + label ? label : row->label, + data_group ? data_group : (row->data_group ? row->data_group : (label ? label : row->label)), + dvodb ? dvodb : NULL, + note ? note : NULL, + image_only, + minidvodb, + minidvodb_group, + minidvodb_name + )) { + if (!psDBRollback(config->dbh)) { + psError(PS_ERR_UNKNOWN, false, "database error sfg"); + } + psError(PS_ERR_UNKNOWN, false, + "failed to trying to queue stage %s %" PRId64,stage, stage_id); + psFree(row); + psFree(output); + return false; + } + psFree(row); + } + } + + + @@ -695,4 +909,7 @@ PXOPT_COPY_S64(config->args, where, "-sky_id", "staticskyRun.sky_id", "=="); PXOPT_COPY_S64(config->args, where, "-skycal_id", "skycalRun.skycal_id", "=="); + PXOPT_COPY_S64(config->args, where, "-diff_id", "diffRun.diff_id", "=="); + PXOPT_COPY_S64(config->args, where, "-ff_id", "fullForceRun.ff_id", "=="); + PXOPT_LOOKUP_STR(stage, config->args, "-stage", false, false); pxcamGetSearchArgs (config, where); // most search arguments based on camera @@ -709,14 +926,20 @@ psString query = NULL; if (strcmp(stage, "cam")==0) { - query = psStringCopy("UPDATE addRun JOIN camRun on cam_id = stage_id JOIN chipRun USING(chip_id) JOIN rawExp USING(exp_id)"); + query = psStringCopy("UPDATE addRun JOIN camRun on (cam_id = stage_id and stage = 'cam') JOIN chipRun USING(chip_id) JOIN rawExp USING(exp_id)"); } if (strcmp(stage, "stack")==0) { - query = psStringCopy("UPDATE addRun JOIN stackRun on stack_id = stage_id"); + query = psStringCopy("UPDATE addRun JOIN stackRun on (stack_id = stage_id and stage = 'stack')"); } if (strcmp(stage, "staticsky")==0) { - query = psStringCopy("UPDATE addRun JOIN staticskyRun on sky_id = stage_id"); + query = psStringCopy("UPDATE addRun JOIN staticskyRun on (sky_id = stage_id and stage = 'staticsky')"); } if (strcmp(stage, "skycal")==0) { - query = psStringCopy("UPDATE addRun JOIN skycalRun on skycal_id = stage_id"); + query = psStringCopy("UPDATE addRun JOIN skycalRun on (skycal_id = stage_id and stage = 'skycal')"); + } + if (strcmp(stage, "diff")==0) { + query = psStringCopy("UPDATE addRun JOIN diffInputSkyfile on (diff_id = stage_id and diff_skyfile_id = stage_extra1 and stage = 'diff') JOIN diffSkyfile on (diff_id, skycell_id) JOIN diffRun using (diff_id)"); + } + if (strcmp(stage, "fullforce")==0) { + query = psStringCopy("UPDATE addRun JOIN fullForceResult on (ff_id = stage_id and warp_id = stage_extra1 and stage = 'fullforce') JOIN fullForceRun on (ff_id)"); } @@ -743,4 +966,6 @@ PXOPT_COPY_S64(config->args, where, "-sky_id", "staticskyRun.sky_id", "=="); PXOPT_COPY_S64(config->args, where, "-skycal_id", "skycalRun.skycal_id", "=="); + PXOPT_COPY_S64(config->args, where, "-ff_id", "fullForceRun.ff_id", "=="); + PXOPT_COPY_S64(config->args, where, "-diff_id", "diffRun.diff_id", "=="); PXOPT_LOOKUP_STR(stage, config->args, "-stage", false, false); pxcamGetSearchArgs (config, where); @@ -764,4 +989,10 @@ query = pxDataGet("addtool_find_pendingexp_skycal.sql"); } + if (strcmp(stage, "diff")==0) { + query = pxDataGet("addtool_find_pendingexp_diff.sql"); + } + if (strcmp(stage, "fullforce")==0) { + query = pxDataGet("addtool_find_pendingexp_ff.sql"); + } @@ -790,4 +1021,15 @@ psStringAppend(&query, " GROUP BY %s", "sky_id, stage_extra1"); } + if (strcmp(stage, "diff") == 0) { + //this group by is needed to join against all the warps (to get camera) + psStringAppend(&query, " GROUP BY %s", "diff_id, stage_extra1"); + } + if (strcmp(stage, "fullforce") == 0) { + //this group by is needed to join against all the warps (to get camera) + psStringAppend(&query, " GROUP BY %s", "ff_id, stage_extra1"); + } + + + // treat limit == 0 as "no limit" if (limit) { @@ -994,4 +1236,9 @@ } else if (strcmp (stage,"skycal") == 0) { query = pxDataGet("addtool_find_processedexp_skycal.sql"); + } else if (strcmp (stage,"diff") == 0) { + query = pxDataGet("addtool_find_processedexp_diff.sql"); + } else if (strcmp (stage,"fullforce") == 0) { + query = pxDataGet("addtool_find_processedexp_ff.sql"); + } else { psFree(where); @@ -1112,4 +1359,10 @@ if (strcmp(stage, "skycal") == 0) { query = pxDataGet("addtool_revertprocessedexp_skycal.sql"); + } + if (strcmp(stage, "diff") == 0) { + query = pxDataGet("addtool_revertprocessedexp_diff.sql"); + } + if (strcmp(stage, "fullforce") == 0) { + query = pxDataGet("addtool_revertprocessedexp_ff.sql"); } Index: branches/eam_branches/ipp-20140904/ippTools/src/addtoolConfig.c =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/src/addtoolConfig.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTools/src/addtoolConfig.c (revision 37804) @@ -55,4 +55,6 @@ psMetadataAddS64(definebyqueryArgs, PS_LIST_TAIL, "-skycal_id", 0, "search by skycal_id", 0); pxcamSetSearchArgs(definebyqueryArgs); + psMetadataAddS64(definebyqueryArgs, PS_LIST_TAIL, "-diff_id", 0, "search by diff_id", 0); + psMetadataAddS64(definebyqueryArgs, PS_LIST_TAIL, "-ff_id", 0, "search by ff_id", 0); psMetadataAddStr(definebyqueryArgs, PS_LIST_TAIL, "-label", PS_META_DUPLICATE_OK, "search by camRun label", NULL); psMetadataAddStr(definebyqueryArgs, PS_LIST_TAIL, "-data_group", PS_META_DUPLICATE_OK, "search by camRun data_group", NULL); @@ -83,4 +85,6 @@ psMetadataAddS64(updaterunArgs, PS_LIST_TAIL, "-sky_id", 0, "search by sky_id", 0); psMetadataAddS64(updaterunArgs, PS_LIST_TAIL, "-skycal_id", 0, "search by skycal_id", 0); + psMetadataAddS64(updaterunArgs, PS_LIST_TAIL, "-diff_id", 0, "search by diff_id",0); + psMetadataAddS64(updaterunArgs, PS_LIST_TAIL, "-ff_id", 0, "search by _ff_id",0); pxcamSetSearchArgs(updaterunArgs); psMetadataAddStr(updaterunArgs, PS_LIST_TAIL, "-label", 0, "search by addRun label", NULL); @@ -101,7 +105,9 @@ psMetadataAddS64(pendingexpArgs, PS_LIST_TAIL, "-add_id", 0, "search by add_id", 0); psMetadataAddS64(pendingexpArgs, PS_LIST_TAIL, "-cam_id", 0, "search by cam_id", 0); - psMetadataAddS64(pendingexpArgs, PS_LIST_TAIL, "-stack_id", 0, "search by stack_id", 0); - psMetadataAddS64(pendingexpArgs, PS_LIST_TAIL, "-sky_id", 0, "search by sky_id", 0); - psMetadataAddS64(pendingexpArgs, PS_LIST_TAIL, "-skycal_id", 0, "search by skycal_id", 0); + psMetadataAddS64(pendingexpArgs, PS_LIST_TAIL, "-stack_id", 0, "search by stack_id", 0); + psMetadataAddS64(pendingexpArgs, PS_LIST_TAIL, "-sky_id", 0, "search by sky_id", 0); + psMetadataAddS64(pendingexpArgs, PS_LIST_TAIL, "-skycal_id", 0, "search by skycal_id", 0); + psMetadataAddS64(pendingexpArgs, PS_LIST_TAIL, "-diff_id", 0, "search by diff_id",0); + psMetadataAddS64(pendingexpArgs, PS_LIST_TAIL, "-ff_id", 0, "search by ff_id",0); pxcamSetSearchArgs(pendingexpArgs); psMetadataAddStr(pendingexpArgs, PS_LIST_TAIL, "-label", PS_META_DUPLICATE_OK, "search by addRun label", NULL); @@ -236,5 +242,5 @@ psMetadata *modes = psMetadataAlloc(); - PXOPT_ADD_MODE("-definebyquery", "create runs from cam stage", ADDTOOL_MODE_DEFINEBYQUERY, definebyqueryArgs); + PXOPT_ADD_MODE("-definebyquery", "create runs from cam/skycal/forcedwarp/diff stage", ADDTOOL_MODE_DEFINEBYQUERY, definebyqueryArgs); PXOPT_ADD_MODE("-updaterun", "change add run properties", ADDTOOL_MODE_UPDATERUN, updaterunArgs); PXOPT_ADD_MODE("-pendingexp", "show pending exps", ADDTOOL_MODE_PENDINGEXP, pendingexpArgs); Index: branches/eam_branches/ipp-20140904/ippTools/src/chiptool.c =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/src/chiptool.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTools/src/chiptool.c (revision 37804) @@ -1454,6 +1454,6 @@ PXOPT_COPY_S64(config->args, where, "-chip_id", "chipRun.chip_id", "=="); PXOPT_COPY_STR(config->args, where, "-reduction", "chipRun.reduction", "=="); - pxAddLabelSearchArgs (config, where, "-label", "label", "=="); - PXOPT_COPY_STR(config->args, where, "-state", "state", "=="); + pxAddLabelSearchArgs (config, where, "-label", "chipRun.label", "=="); + PXOPT_COPY_STR(config->args, where, "-state", "chipRun.state", "=="); psString query = pxDataGet("chiptool_run.sql"); Index: branches/eam_branches/ipp-20140904/ippTools/src/fftool.c =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/src/fftool.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTools/src/fftool.c (revision 37804) @@ -859,5 +859,5 @@ PXOPT_COPY_STR(config->args, where, "-skycell_id", "stackRun.skycell_id", "LIKE"); PXOPT_COPY_STR(config->args, where, "-filter", "stackRun.filter", "LIKE"); - PXOPT_COPY_S16(config->args, where, "-fault", "staticskyResult.fault", "=="); + PXOPT_COPY_S16(config->args, where, "-fault", "fullForceSummary.fault", "=="); pxskycellAddWhere(config, where); Index: branches/eam_branches/ipp-20140904/ippTools/src/laptool.c =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/src/laptool.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTools/src/laptool.c (revision 37804) @@ -25,4 +25,5 @@ static bool pendingrunMode(pxConfig *config); static bool updaterunMode(pxConfig *config); +static bool revertsassMode(pxConfig *config); // Exposure level @@ -69,5 +70,6 @@ MODECASE(LAPTOOL_MODE_PENDINGRUN, pendingrunMode); MODECASE(LAPTOOL_MODE_UPDATERUN, updaterunMode); - + MODECASE(LAPTOOL_MODE_REVERTSASS, revertsassMode); + MODECASE(LAPTOOL_MODE_PENDINGEXP, pendingexpMode); MODECASE(LAPTOOL_MODE_EXPOSURES, exposuresMode); @@ -695,4 +697,51 @@ } +static bool revertsassMode(pxConfig *config) +{ + PS_ASSERT_PTR_NON_NULL(config, false); + + psMetadata *where = psMetadataAlloc(); + PXOPT_COPY_S64(config->args, where, "-lap_id", "lap_id", "=="); + PXOPT_COPY_S64(config->args, where, "-sass_id", "final_sass_id", "=="); + if (!psListLength(where->list)) { + psFree(where); + psError(PXTOOLS_ERR_CONFIG, false, "search parameters are required"); + return false; + } + + psString query = pxDataGet("laptool_revertsass.sql"); + if (!query) { + psError(PXTOOLS_ERR_SYS, false, "failed to retrieve SQL statement"); + return(false); + } + psString whereClause = psDBGenerateWhereConditionSQL(where, NULL); + psStringAppend(&query, "\n WHERE %s", whereClause); + psFree(whereClause); + + if (!psDBTransaction(config->dbh)) { + psError(PS_ERR_UNKNOWN, false, "database error"); + return false; + } + + if (!p_psDBRunQuery(config->dbh, query)) { + psError(PS_ERR_UNKNOWN, false, "database error"); + psFree(query); + return false; + } + psFree(query); + + long numUpdated = psDBAffectedRows(config->dbh); + if (numUpdated != 1) { + psError(PS_ERR_UNKNOWN, false, "should have updated only one row. updated %ld", numUpdated); + return false; + } + + if (!psDBCommit(config->dbh)) { + psError(PS_ERR_UNKNOWN, false, "database error"); + return false; + } + + return true; +} // Exposure level Index: branches/eam_branches/ipp-20140904/ippTools/src/laptool.h =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/src/laptool.h (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTools/src/laptool.h (revision 37804) @@ -19,4 +19,5 @@ LAPTOOL_MODE_EXPOSURES, LAPTOOL_MODE_STACKS, + LAPTOOL_MODE_REVERTSASS, LAPTOOL_MODE_UPDATEEXP, LAPTOOL_MODE_DIFFCHECK, Index: branches/eam_branches/ipp-20140904/ippTools/src/laptoolConfig.c =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/src/laptoolConfig.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTools/src/laptoolConfig.c (revision 37804) @@ -117,4 +117,9 @@ ADD_OPT(S64, updaterunArgs, "-set_final_sass_id", "set final stack sass_id", 0); + // -revertsass + psMetadata *revertsassArgs = psMetadataAlloc(); + ADD_OPT(S64, revertsassArgs, "-lap_id", "search by lap run id", 0); + ADD_OPT(S64, revertsassArgs, "-sass_id", "sass id that is being removed", 0); + // -pendingexp psMetadata *pendingexpArgs = psMetadataAlloc(); @@ -136,5 +141,5 @@ ADD_OPT(Bool,stacksArgs, "-simple", "use the simple output format", false); ADD_OPT(U64, stacksArgs, "-limit", "limit result set to N items", 0); - + // -updateexp psMetadata *updateexpArgs = psMetadataAlloc(); @@ -227,4 +232,5 @@ PXOPT_ADD_MODE("-listrun", "", LAPTOOL_MODE_PENDINGRUN, listrunArgs); PXOPT_ADD_MODE("-updaterun", "", LAPTOOL_MODE_UPDATERUN, updaterunArgs); + PXOPT_ADD_MODE("-revertsass", "mode to remove the sass_id", LAPTOOL_MODE_REVERTSASS, revertsassArgs); PXOPT_ADD_MODE("-pendingexp", "", LAPTOOL_MODE_PENDINGEXP, pendingexpArgs); PXOPT_ADD_MODE("-exposures", "", LAPTOOL_MODE_EXPOSURES, exposuresArgs); Index: branches/eam_branches/ipp-20140904/ippTools/src/pxadd.c =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/src/pxadd.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTools/src/pxadd.c (revision 37804) @@ -164,4 +164,14 @@ psMemSetPersistent(query, true); } + else if (strcmp(stage,"diff") == 0) { + query = pxDataGet("addtool_queue_diff_id.sql"); + psMemSetPersistent(query, true); + } + else if (strcmp(stage,"fullforce") == 0) { + query = pxDataGet("addtool_queue_ff_id.sql"); + psMemSetPersistent(query, true); + } + + else { Index: branches/eam_branches/ipp-20140904/ippTools/src/stacktool.c =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/src/stacktool.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTools/src/stacktool.c (revision 37804) @@ -41,4 +41,5 @@ static bool sumskyfileMode(pxConfig *config); static bool sassskyfileMode(pxConfig *config); +static bool updatesassMode(pxConfig *config); static bool revertsumskyfileMode(pxConfig *config); static bool tosummaryMode(pxConfig *config); @@ -82,5 +83,6 @@ MODECASE(STACKTOOL_MODE_ADDSUMSKYFILE, addsumskyfileMode); MODECASE(STACKTOOL_MODE_SUMSKYFILE, sumskyfileMode); - MODECASE(STACKTOOL_MODE_SASSSKYFILE, sassskyfileMode); + MODECASE(STACKTOOL_MODE_SASSSKYFILE, sassskyfileMode); + MODECASE(STACKTOOL_MODE_UPDATESASS, updatesassMode); MODECASE(STACKTOOL_MODE_REVERTSUMSKYFILE, revertsumskyfileMode); MODECASE(STACKTOOL_MODE_TOSUMMARY, tosummaryMode); @@ -756,6 +758,12 @@ PXOPT_COPY_STR(config->args, where, "-data_group", "stackRun.data_group", "=="); PXOPT_COPY_STR(config->args, where, "-state", "stackRun.state", "=="); - PXOPT_COPY_S64(config->args, where, "-sass_id", "stackAssociationMap.sass_id", "=="); + PXOPT_LOOKUP_S16(fault, config->args, "-fault", false, false); + PXOPT_LOOKUP_S64(sass_id, config->args, "-sass_id", false, false); + + if (sass_id) { + PXOPT_COPY_S64(config->args, where, "-sass_id", "stackAssociationMap.sass_id", "=="); + } + if (!psListLength(where->list)) { psFree(where); @@ -771,5 +779,8 @@ PXOPT_COPY_S16(config->args, where, "-fault", "stackSumSkyfile.fault", "=="); } - + if (sass_id) { + psStringAppend(&query, " JOIN stackAssociationMap USING(stack_id)"); + } + // pxUpdateRun gets parameters from config->args and updates bool result = pxUpdateRun(config, where, &query, "stackRun", "stack_id", "stackSumSkyfile", true, false); @@ -1295,9 +1306,4 @@ PXOPT_COPY_STR(config->args, where, "-data_group", "stackAssociation.data_group", "LIKE"); -// The following three selectors are incompatible with the sql so omit them -// PXOPT_COPY_S64(config->args, where, "-warp_id", "warpRun.warp_id", "=="); -// PXOPT_COPY_S64(config->args, where, "-exp_id", "rawExp.exp_id", "=="); -// PXOPT_COPY_STR(config->args, where, "-exp_name", "rawExp.exp_name", "=="); - PXOPT_LOOKUP_BOOL(all, config->args, "-all", false); @@ -1369,4 +1375,49 @@ } +static bool updatesassMode(pxConfig *config) +{ + PS_ASSERT_PTR_NON_NULL(config, false); + + psMetadata *where = psMetadataAlloc(); + PXOPT_COPY_S64(config->args, where, "-sass_id", "sass_id", "=="); + PXOPT_COPY_STR(config->args, where, "-data_group", "data_group", "=="); + PXOPT_COPY_STR(config->args, where, "-projection_cell", "projection_cell", "=="); + PXOPT_COPY_STR(config->args, where, "-tess_id", "tess_id", "=="); + PXOPT_COPY_STR(config->args, where, "-filter", "filter", "=="); + if (!psListLength(where->list)) { + psFree(where); + psError(PXTOOLS_ERR_CONFIG, false, "search parameters are required"); + return false; + } + + psMetadata *values = psMetadataAlloc(); + PXOPT_COPY_STR(config->args, values, "-set_data_group", "data_group", "=="); + if (!psListLength(values->list)) { + psFree(values); + psError(PXTOOLS_ERR_CONFIG, false, "data_group to assign not found"); + return false; + } + + if (!psDBTransaction(config->dbh)) { + psError(PS_ERR_UNKNOWN, false, "database error"); + return false; + } + + long rows = psDBUpdateRows(config->dbh, "stackAssociation", where, values); + psFree(values); + psFree(where); + + if (rows != 1) { + psError(PS_ERR_UNKNOWN, false, "should have updated only one row. updated %ld", rows); + return false; + } + + if (!psDBCommit(config->dbh)) { + psError(PS_ERR_UNKNOWN, false, "database error"); + return false; + } + + return true; +} static bool revertsumskyfileMode(pxConfig *config) Index: branches/eam_branches/ipp-20140904/ippTools/src/stacktool.h =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/src/stacktool.h (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTools/src/stacktool.h (revision 37804) @@ -35,4 +35,5 @@ STACKTOOL_MODE_SUMSKYFILE, STACKTOOL_MODE_SASSSKYFILE, + STACKTOOL_MODE_UPDATESASS, STACKTOOL_MODE_REVERTSUMSKYFILE, STACKTOOL_MODE_TOSUMMARY, Index: branches/eam_branches/ipp-20140904/ippTools/src/stacktoolConfig.c =================================================================== --- branches/eam_branches/ipp-20140904/ippTools/src/stacktoolConfig.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ippTools/src/stacktoolConfig.c (revision 37804) @@ -238,4 +238,14 @@ psMetadataAddBool(sassskyfileArgs, PS_LIST_TAIL, "-all", 0, "enable search without arguments", false); + // -updatesass + psMetadata *updatesassArgs = psMetadataAlloc(); + psMetadataAddS64(updatesassArgs, PS_LIST_TAIL, "-sass_id", 0, "search by stack association id", 0); + psMetadataAddStr(updatesassArgs, PS_LIST_TAIL, "-data_group", 0, "search by data group", 0); + psMetadataAddStr(updatesassArgs, PS_LIST_TAIL, "-projection_cell", 0, "search by projection cell", 0); + psMetadataAddStr(updatesassArgs, PS_LIST_TAIL, "-tess_id", 0, "search by tess_id", 0); + psMetadataAddStr(updatesassArgs, PS_LIST_TAIL, "-filter", 0, "search by filter", 0); + + psMetadataAddStr(updatesassArgs, PS_LIST_TAIL, "-set_data_group", 0, "data_group to assign", 0); + // -revertsumskyfile psMetadata *revertsumskyfileArgs= psMetadataAlloc(); @@ -348,5 +358,6 @@ PXOPT_ADD_MODE("-donecleanup", "show runs that have been cleaned", STACKTOOL_MODE_DONECLEANUP, donecleanupArgs); - PXOPT_ADD_MODE("-sassskyfile", "list results of stackAssociation", STACKTOOL_MODE_SASSSKYFILE, sassskyfileArgs); + PXOPT_ADD_MODE("-sassskyfile", "list results of stackAssociation", STACKTOOL_MODE_SASSSKYFILE, sassskyfileArgs); + PXOPT_ADD_MODE("-updatesass", "update data_group of stackAssociation", STACKTOOL_MODE_UPDATESASS, updatesassArgs); PXOPT_ADD_MODE("-tosummary", "show runs that can be summarized", STACKTOOL_MODE_TOSUMMARY, tosummaryArgs); PXOPT_ADD_MODE("-summary", "show runs that have been summarized", STACKTOOL_MODE_SUMMARY, summaryArgs); Index: branches/eam_branches/ipp-20140904/mops/mops_ringsv3/populate_ringsv3_database.py =================================================================== --- branches/eam_branches/ipp-20140904/mops/mops_ringsv3/populate_ringsv3_database.py (revision 37804) +++ branches/eam_branches/ipp-20140904/mops/mops_ringsv3/populate_ringsv3_database.py (revision 37804) @@ -0,0 +1,82 @@ +#!/usr/bin/env python2.7 + +import MySQLdb +import pyfits +import numpy + +################################################################### +# "Parameters" +################################################################### +# Name of the table which is going to be created and populated +tableName = "skycells" +# Name of the database host +dbhost = 'localhost' +dbuser = 'mops' +dbpasswd = 'mops' +dbname = 'ringsV3' +# Name of the tesselation file +tesselation = '/local/ipp/gpc1/tess/RINGS.V3/Images.dat' +################################################################### + +def create_database(db): + cursor = db.cursor() + statement = "DROP TABLE IF EXISTS %s;" % tableName + cursor.execute(statement) + + statement = """CREATE TABLE %s ( + skycell CHAR(16) NOT NULL PRIMARY KEY COMMENT 'Skycell ID (aka Name in Images.dat file)', + raProjectionCenter DOUBLE NOT NULL COMMENT 'Right Ascension of the center of projection cell (in degrees, aka crval1)', + craProjectionCenter DOUBLE NOT NULL COMMENT 'Cosine of ra of center of projection cell', + sraProjectionCenter DOUBLE NOT NULL COMMENT 'Sine of ra of center of projection cell', + decProjectionCenter DOUBLE NOT NULL COMMENT 'Declination of center of projection cell (in degrees, aka crval2)', + cdecProjectionCenter DOUBLE NOT NULL COMMENT 'Cosine of decl of center of projection cell', + sdecProjectionCenter DOUBLE NOT NULL COMMENT 'Sine of decl of center of projection cell', + xSkycellCenter DOUBLE NOT NULL COMMENT 'X coordinate of the skycell center (in pixels, aka crpix1)', + ySkycellCenter DOUBLE NOT NULL COMMENT 'Y coordinate of the skycell center (in pixels, aka crpix2)' +) ENGINE=InnoDB COMMENT='Contents of tesselation used for MOPS diffs mapped to a database table. Notes: (a) Scales on both X and Y axis (aka cdelt1, cdelt2) are 5 arcseconds. (b) Rotation coordinates (aka pc[12]_[12]) in from (ra,dec) coordinates to image plane coordinates are ((-1,0),(0,1)). (c) The database has been populated on 2014-10-20 by schastel using /local/ipp/gpc1/tess/RINGS.V3/Images.dat whose md5sum is 8eb517fe2e894e66488d154f5bb7c030; (d) Associated software: SVN repository: https://svn.pan-starrs.ifa.hawaii.edu/repo/ipp/trunk/mops/mops_ringsv3' ; +""" % (tableName) + cursor.execute(statement) + cursor.close() + +def populate_database(database): + fits = pyfits.open(tesselation) + table = fits[1].data + raProjectionCenterRadians = numpy.radians(table.field('CRVAL1')) + craProjectionCenter = numpy.cos(raProjectionCenterRadians) + sraProjectionCenter = numpy.sin(raProjectionCenterRadians) + decProjectionCenterRadians = numpy.radians(table.field('CRVAL2')) + cdecProjectionCenter = numpy.cos(decProjectionCenterRadians) + sdecProjectionCenter = numpy.sin(decProjectionCenterRadians) + data = [ table.field('NAME'), + table.field('CRVAL1'), craProjectionCenter, sraProjectionCenter, + table.field('CRVAL2'), cdecProjectionCenter, sdecProjectionCenter, + table.field('CRPIX1'), table.field('CRPIX2'), + ] + query = """INSERT INTO %s (skycell, +raProjectionCenter, craProjectionCenter, sraProjectionCenter, +decProjectionCenter, cdecProjectionCenter, sdecProjectionCenter, +xSkycellCenter, ySkycellCenter) VALUES (%%s, %%s, %%s, %%s, %%s, %%s, %%s, %%s, %%s)""" % tableName + # print query + remapped = zip(*data) + # print + # print query % remapped[0] + # print + # print remapped[0:5] + begin = 0 + step = 10000 + while begin < len(remapped): + cursor = database.cursor() + cursor.executemany(query, remapped[begin:begin+step]) + cursor.close() + database.commit() + begin += step + cursor.executemany(query, remapped[begin:-1]) + fits.close() + cursor.close() + +if __name__ == "__main__": + database = MySQLdb.connect(dbhost, dbuser, dbpasswd, dbname) + create_database(database) + populate_database(database) + database.commit() + database.close() Index: branches/eam_branches/ipp-20140904/ppImage/src/ppImageAddNoise.c =================================================================== --- branches/eam_branches/ipp-20140904/ppImage/src/ppImageAddNoise.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ppImage/src/ppImageAddNoise.c (revision 37804) @@ -18,4 +18,6 @@ return true; } + + // psTimerStart("add.noise"); // grizy variances to add to turn MD exposure -> 3pi, calculated from the DRM @@ -96,4 +98,6 @@ psFree(stats); + // psLogMsg ("ppImage", 5, "add noise: %f sec\n", psTimerMark ("add.noise")); + return true; } Index: branches/eam_branches/ipp-20140904/ppImage/src/ppImageAuxiliaryMask.c =================================================================== --- branches/eam_branches/ipp-20140904/ppImage/src/ppImageAuxiliaryMask.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ppImage/src/ppImageAuxiliaryMask.c (revision 37804) @@ -67,4 +67,6 @@ return true; } + + // psTimerStart("aux.mask"); bool status; @@ -208,4 +210,6 @@ } + // psLogMsg ("ppImage", 5, "auxiliary mask: %f sec\n", psTimerMark ("aux.mask")); + return true; } Index: branches/eam_branches/ipp-20140904/ppImage/src/ppImageCheckCTE.c =================================================================== --- branches/eam_branches/ipp-20140904/ppImage/src/ppImageCheckCTE.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ppImage/src/ppImageCheckCTE.c (revision 37804) @@ -17,4 +17,6 @@ return true; } + + // psTimerStart("check.cte"); // add recipe options supplied on command line @@ -106,4 +108,6 @@ psFree (statsDefaults); + // psLogMsg ("ppImage", 5, "check CTE: %f sec\n", psTimerMark ("check.cte")); + return true; } Index: branches/eam_branches/ipp-20140904/ppImage/src/ppImageCheckNoise.c =================================================================== --- branches/eam_branches/ipp-20140904/ppImage/src/ppImageCheckNoise.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ppImage/src/ppImageCheckNoise.c (revision 37804) @@ -17,4 +17,6 @@ return true; } + + // psTimerStart("check.noise"); // add recipe options supplied on command line @@ -124,4 +126,6 @@ psFree (statsDefaults); + // psLogMsg ("ppImage", 5, "check noise: %f sec\n", psTimerMark ("check.noise")); + return true; } Index: branches/eam_branches/ipp-20140904/ppImage/src/ppImageDetrendFringe.c =================================================================== --- branches/eam_branches/ipp-20140904/ppImage/src/ppImageDetrendFringe.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ppImage/src/ppImageDetrendFringe.c (revision 37804) @@ -327,4 +327,6 @@ pmCell *cell = NULL; + // psTimerStart("apply.fringe"); + assert (options->doFringe); // do not call if not needed assert (inputView->chip != -1); @@ -386,4 +388,6 @@ } + // psLogMsg ("ppImage", 5, "apply fringe: %f sec\n", psTimerMark ("apply.fringe")); + psFree (view); return true; Index: branches/eam_branches/ipp-20140904/ppImage/src/ppImageDetrendPattern.c =================================================================== --- branches/eam_branches/ipp-20140904/ppImage/src/ppImageDetrendPattern.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ppImage/src/ppImageDetrendPattern.c (revision 37804) @@ -35,5 +35,7 @@ // chip. + // psTimerStart("apply.pattern"); if (options->doPatternRow) { + bool status; pmHDU *hdu = pmHDUFromChip(chip); @@ -193,4 +195,7 @@ pattern_done: + + // psLogMsg ("ppImage", 5, "apply pattern: %f sec\n", psTimerMark ("apply.pattern")); + return(true); } Index: branches/eam_branches/ipp-20140904/ppImage/src/ppImageDetrendReadout.c =================================================================== --- branches/eam_branches/ipp-20140904/ppImage/src/ppImageDetrendReadout.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ppImage/src/ppImageDetrendReadout.c (revision 37804) @@ -7,4 +7,6 @@ bool ppImageDetrendReadout(pmConfig *config, ppImageOptions *options, pmFPAview *view) { + // psTimerStart("detrend.readout"); + // construct a view for the detrend images (which have only one readout) pmFPAview *detview = pmFPAviewAlloc(0); @@ -35,4 +37,5 @@ } } + // psLogMsg ("ppImage", 6, "check gain: %f sec\n", psTimerMark ("detrend.readout")); // Check to see if we're in a chip that contains video @@ -56,4 +59,5 @@ } done_video_check: + // psLogMsg ("ppImage", 6, "check video: %f sec\n", psTimerMark ("detrend.readout")); // Masking on the basis of pixel value needs to be done before anything else, so the values are pristine. @@ -66,4 +70,5 @@ return false; } + // psLogMsg ("ppImage", 6, "generate mask: %f sec\n", psTimerMark ("detrend.readout")); } // apply the externally supplied mask to the input->mask pixels @@ -77,4 +82,5 @@ } pmMaskBadPixels(input, mask, options->maskValue); + // psLogMsg ("ppImage", 6, "apply mask: %f sec\n", psTimerMark ("detrend.readout")); } @@ -83,4 +89,5 @@ // extern bool ppImageBurntoolApply(pmConfig *, ppImageOptions *, pmFPAview *, pmReadout *); ppImageBurntoolApply(config, options, view, input); + // psLogMsg ("ppImage", 6, "apply burntool: %f sec\n", psTimerMark ("detrend.readout")); } @@ -93,4 +100,5 @@ ppImageBurntoolMask(config,options,view,input); } + // psLogMsg ("ppImage", 6, "apply burntool mask: %f sec\n", psTimerMark ("detrend.readout")); } @@ -102,4 +110,5 @@ return false; } + // psLogMsg ("ppImage", 6, "subtract overscan: %f sec\n", psTimerMark ("detrend.readout")); } @@ -111,4 +120,5 @@ return(false); } + // psLogMsg ("ppImage", 6, "nonlinear correction: %f sec\n", psTimerMark ("detrend.readout")); } @@ -148,4 +158,5 @@ return false; } + // psLogMsg ("ppImage", 6, "apply bias: %f sec\n", psTimerMark ("detrend.readout")); } @@ -218,4 +229,5 @@ } psFree (noiseImage); + // psLogMsg ("ppImage", 6, "generate variance: %f sec\n", psTimerMark ("detrend.readout")); } @@ -226,4 +238,5 @@ return false; } + // psLogMsg ("ppImage", 6, "apply dark: %f sec\n", psTimerMark ("detrend.readout")); } @@ -235,4 +248,5 @@ return false; } + // psLogMsg ("ppImage", 6, "mask remnance: %f sec\n", psTimerMark ("detrend.readout")); } @@ -244,4 +258,5 @@ return false; } + // psLogMsg ("ppImage", 6, "shutter correction: %f sec\n", psTimerMark ("detrend.readout")); } @@ -253,4 +268,5 @@ return false; } + // psLogMsg ("ppImage", 6, "apply flat: %f sec\n", psTimerMark ("detrend.readout")); } @@ -276,4 +292,5 @@ psBinaryOp(input->image, input->image, "*", psScalarAlloc(norm, PS_TYPE_F32)); + // psLogMsg ("ppImage", 6, "renormalize: %f sec\n", psTimerMark ("detrend.readout")); } @@ -324,4 +341,6 @@ psFree (classID); + + // psLogMsg ("ppImage", 6, "renormalize by class: %f sec\n", psTimerMark ("detrend.readout")); } # endif @@ -333,8 +352,10 @@ return false; } + // psLogMsg ("ppImage", 6, "measure fringe: %f sec\n", psTimerMark ("detrend.readout")); } ppImageMemoryDump("detrend"); + // psLogMsg ("ppImage", 5, "detrend readout: %f sec\n", psTimerMark ("detrend.readout")); psFree(detview); return true; Index: branches/eam_branches/ipp-20140904/ppImage/src/ppImageMosaic.c =================================================================== --- branches/eam_branches/ipp-20140904/ppImage/src/ppImageMosaic.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ppImage/src/ppImageMosaic.c (revision 37804) @@ -9,4 +9,6 @@ { bool status; // Status of MD lookup + + // psTimerStart("mosaic.chip"); pmFPAfile *in = psMetadataLookupPtr(&status, config->files, inFile); // Input file @@ -63,4 +65,6 @@ } + // psLogMsg ("ppImage", 5, "mosaic chip: %f sec\n", psTimerMark ("mosaic.chip")); + return status; } Index: branches/eam_branches/ipp-20140904/ppImage/src/ppImageParityFlip.c =================================================================== --- branches/eam_branches/ipp-20140904/ppImage/src/ppImageParityFlip.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ppImage/src/ppImageParityFlip.c (revision 37804) @@ -19,4 +19,6 @@ if (!options->applyParity) return true; + + // psTimerStart("parity.flip"); // find the currently selected readout @@ -152,4 +154,5 @@ } } + // psLogMsg ("ppImage", 5, "parity flip: %f sec\n", psTimerMark ("parity.flip")); FIX_CONCEPT(cell->concepts, "CELL.XPARITY", S32, xParityTarget); Index: branches/eam_branches/ipp-20140904/ppImage/src/ppImagePixelStats.c =================================================================== --- branches/eam_branches/ipp-20140904/ppImage/src/ppImagePixelStats.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ppImage/src/ppImagePixelStats.c (revision 37804) @@ -10,4 +10,6 @@ { bool mdok; // Status of MD lookup + + // psTimerStart("pixel.stats"); // loop over the cells/readouts for this chip @@ -81,4 +83,6 @@ } + // psLogMsg ("ppImage", 5, "measure pixel stats: %f sec\n", psTimerMark ("pixel.stats")); + psFree (view); return true; Index: branches/eam_branches/ipp-20140904/ppImage/src/ppImageRebinReadout.c =================================================================== --- branches/eam_branches/ipp-20140904/ppImage/src/ppImageRebinReadout.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ppImage/src/ppImageRebinReadout.c (revision 37804) @@ -12,4 +12,6 @@ pmFPAfile *outFile = psMetadataLookupPtr(NULL, config->files, outName); if (outFile == NULL) return false; + + // psTimerStart("rebin.chip"); // XXX double check that chip != -1? @@ -40,4 +42,6 @@ } + // psLogMsg ("ppImage", 5, "rebin chip for %s: %f sec\n", outName, psTimerMark ("rebin.chip")); + return true; } Index: branches/eam_branches/ipp-20140904/ppImage/src/ppImageSetMaskBits.c =================================================================== --- branches/eam_branches/ipp-20140904/ppImage/src/ppImageSetMaskBits.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ppImage/src/ppImageSetMaskBits.c (revision 37804) @@ -11,4 +11,6 @@ return false; } + + // psTimerStart("set.mask.bits"); // at this point we know we have valid values for required entries SAT, BAD, FLAT, BLANK: @@ -57,4 +59,6 @@ options->markValue); + // psLogMsg ("ppImage", 5, "set mask bits: %f sec\n", psTimerMark ("set.mask.bits")); + return true; } Index: branches/eam_branches/ipp-20140904/ppTranslate/src/ppMops.h =================================================================== --- branches/eam_branches/ipp-20140904/ppTranslate/src/ppMops.h (revision 37803) +++ branches/eam_branches/ipp-20140904/ppTranslate/src/ppMops.h (revision 37804) @@ -6,5 +6,5 @@ #define IN_EXTNAME "SkyChip.psf" // Extension name for data in input -#define OBSERVATORY_CODE "F51" // IAU Observatory Code +// #define OBSERVATORY_CODE "F51" // IAU Observatory Code #define OUT_EXTNAME "MOPS_TRANSIENT_DETECTIONS" // Extension name for data in output #define SOURCE_MASK (PM_SOURCE_MODE_FAIL | PM_SOURCE_MODE_BADPSF | PM_SOURCE_MODE_SATURATED | \ @@ -23,4 +23,6 @@ psS64 warp_id; // Warp stage identifier psS64 diff_id; // Diff stage identifier + psString camera; // Camera name + psString obscode; // observatory code (derived from camera) bool positive; // Sense of subtraction, T=positive, F=negative float zp, zpErr; // Magnitude zero point and error @@ -59,5 +61,6 @@ psVector *raExtErr, *decExtErr; // Error in sky coordinates for RA_EXT and DEC_EXT psVector *mask; // Mask for detections - float platescale; //Plate scale at centroid + float platescale; // Plate scale at centroid + psString fpashutoutc; // FPA shutoutc } ppMopsDetections; Index: branches/eam_branches/ipp-20140904/ppTranslate/src/ppMopsArguments.c =================================================================== --- branches/eam_branches/ipp-20140904/ppTranslate/src/ppMopsArguments.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ppTranslate/src/ppMopsArguments.c (revision 37804) @@ -42,4 +42,6 @@ static void mopsArgumentsFree(ppMopsArguments *args) { + psFree(args->camera); + psFree(args->obscode); psFree(args->input); psFree(args->exp_name); @@ -63,4 +65,6 @@ args->warp_id = 0; args->diff_id = 0; + args->camera = NULL; + args->obscode = NULL; args->zp = NAN; args->positive = true; @@ -94,4 +98,5 @@ psMetadataAddS64(arguments, PS_LIST_TAIL, "-warp_id", 0, "Warp stage identifier", 0); psMetadataAddS64(arguments, PS_LIST_TAIL, "-diff_id", 0, "Diff stage identifier", 0); + psMetadataAddStr(arguments, PS_LIST_TAIL, "-camera", 0, "Camera name", "No comment provided"); psMetadataAddBool(arguments, PS_LIST_TAIL, "-inverse", 0, "Inverse subtraction?", false); psMetadataAddF32(arguments, PS_LIST_TAIL, "-zp", 0, "Magnitude zero point", NAN); @@ -126,4 +131,5 @@ args->warp_id = psMetadataLookupS64(NULL, arguments, "-warp_id"); args->diff_id = psMetadataLookupS64(NULL, arguments, "-diff_id"); + args->camera = psMemIncrRefCounter(psMetadataLookupStr(NULL, arguments, "-camera")); args->positive = !psMetadataLookupBool(NULL, arguments, "-inverse"); // NOTE: negated args->zp = psMetadataLookupF32(NULL, arguments, "-zp"); @@ -137,4 +143,19 @@ args->shutoutc = psMemIncrRefCounter(psMetadataLookupStr(NULL, arguments, "-shutoutc")); + if (!args->camera) { + psError(PS_ERR_BAD_PARAMETER_VALUE, true, "camera name not provided."); + return NULL; + } + if (!strcasecmp(args->camera, "GPC1")) { + args->obscode = psStringCopy("F51"); + } + if (!strcasecmp(args->camera, "GPC2")) { + args->obscode = psStringCopy("F52"); + } + if (!args->obscode) { + psError(PS_ERR_BAD_PARAMETER_VALUE, true, "unknown camera name."); + return NULL; + } + psTrace("ppMops.args", 1, "Done parsing command-line arguments\n"); Index: branches/eam_branches/ipp-20140904/ppTranslate/src/ppMopsDetections.c =================================================================== --- branches/eam_branches/ipp-20140904/ppTranslate/src/ppMopsDetections.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ppTranslate/src/ppMopsDetections.c (revision 37804) @@ -20,4 +20,5 @@ psFree(det->decErr); psFree(det->mask); + psFree(det->fpashutoutc); return; } @@ -48,4 +49,5 @@ det->mask = NULL; det->diffSkyfileId = 0; + det->fpashutoutc = NULL; return det; } Index: branches/eam_branches/ipp-20140904/ppTranslate/src/ppMopsRead.c =================================================================== --- branches/eam_branches/ipp-20140904/ppTranslate/src/ppMopsRead.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ppTranslate/src/ppMopsRead.c (revision 37804) @@ -105,4 +105,5 @@ det->naxis1 = psMetadataLookupS32(NULL, header, "IMNAXIS1"); // Number of columns det->naxis2 = psMetadataLookupS32(NULL, header, "IMNAXIS2"); // Number of rows + det->fpashutoutc = psMemIncrRefCounter(psMetadataLookupStr(NULL, header, "FPA.SHUTOUTC")); psFree(header); @@ -241,7 +242,8 @@ long numGood = 0; // Number of good rows for (long row = 0; row < size; row++) { - psU32 flags = flagsV->data.U32[row]; // psFitsTableGetU32(NULL, table, row, "FLAGS"); + //psU32 flags = flagsV->data.U32[row]; // psFitsTableGetU32(NULL, table, row, "FLAGS"); + psU64 flags = flagsV->data.U64[row]; // table reads in as U64 if (flags & SOURCE_MASK) { - psTrace("ppMops.read", 10, "Discarding row %ld from input %d because of flags: %ud", row, i, flags); + psTrace("ppMops.read", 10, "Discarding row %ld from input %d because of flags: %lud", row, i, flags); det->mask->data.U8[row] = 0xFF; continue; Index: branches/eam_branches/ipp-20140904/ppTranslate/src/ppMopsWrite.c =================================================================== --- branches/eam_branches/ipp-20140904/ppTranslate/src/ppMopsWrite.c (revision 37803) +++ branches/eam_branches/ipp-20140904/ppTranslate/src/ppMopsWrite.c (revision 37804) @@ -89,5 +89,5 @@ psWarning("no inputs with surviving detections. output header will be incomplete"); } - psMetadataAddStr(header, PS_LIST_TAIL, "OBSCODE", 0, "IAU Observatory code", OBSERVATORY_CODE); + psMetadataAddStr(header, PS_LIST_TAIL, "OBSCODE", 0, "IAU Observatory code", args->obscode); psMetadataAddF32(header, PS_LIST_TAIL, "MAGZP", 0, "Magnitude zero point", args->zp); //MOPS want the name ZEROPOINT for MAGZP @@ -102,4 +102,5 @@ psMetadataAddF32(header, PS_LIST_TAIL, "PLTSCALE_EXT", 0, "Plate scale at centroid", det->platescale); psMetadataAddF32(header, PS_LIST_TAIL, "PLTSCALE", 0, "Plate scale at centroid", det->platescale); + psMetadataAddStr(header, PS_LIST_TAIL, "FPA.SHUTOUTC", 0, "Time of exposure", det->fpashutoutc); //field in header that tells about the CMF version Index: branches/eam_branches/ipp-20140904/psLib/src/imageops/psImageMapFit.c =================================================================== --- branches/eam_branches/ipp-20140904/psLib/src/imageops/psImageMapFit.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psLib/src/imageops/psImageMapFit.c (revision 37804) @@ -18,4 +18,5 @@ #include +#include #include "psError.h" #include "psAbort.h" @@ -35,4 +36,5 @@ #include "psImageStructManip.h" #include "psImageMap.h" +#include "psSparse.h" // #include "psImagePixelInterpolate.h" // #include "psImageUnbin.h" @@ -118,5 +120,5 @@ psImageInit (A, 0.0); psVectorInit (B, 0.0); - + // we are looping over the Nx,Ny image map elements; // the matrix equation contains Nx*Ny rows and columns @@ -262,5 +264,5 @@ int I = n + Nx * m; B->data.F32[I] = fi_rx_ry + fi_rx_py + fi_px_ry + fi_px_py; - + // insert these values into their corresponding locations in A, B // float Sum = 0.0; @@ -273,4 +275,5 @@ int J = (n + jn) + Nx * (m + jm); A->data.F32[J][I] = sA[jn][jm]; + // fprintf (stderr, "A %d %d (%d %d : %d %d): %f\n", I, J, n, m, n + jn, m + jm, sA[jn][jm]); // Sum += sA[jn][jm]; @@ -320,5 +323,5 @@ return true; } - + // set bad values to NaN for (int i = 0; i < Nx*Ny; i++) { @@ -341,5 +344,4 @@ psFree (B); psFree (Empty); - *pGoodFit = true; return true; @@ -402,4 +404,432 @@ psS32 Nkeep = 0; if (!psImageMapFit(pGoodFit, map, mask, maskValue, x, y, f, df)) { + psError(PS_ERR_UNKNOWN, false, "Could not fit image map.\n"); + psFree(resid); + if (!inMask) psFree (mask); + return false; + } + if (!*pGoodFit) { + psWarning ("bad fit to image map, try something else"); + psFree(resid); + if (!inMask) psFree (mask); + return true; + } + + psVector *fit = psImageMapEvalVector(map, mask, maskValue, x, y); + if (fit == NULL) { + psError(PS_ERR_UNKNOWN, false, "Failure in psImageMapEvalVector().\n"); + psFree(resid); + if (!inMask) psFree (mask); + return false; + } + for (int i = 0 ; i < f->n ; i++) { + resid->data.F32[i] = (f->data.F32[i] - fit->data.F32[i]); + } + + if (!psVectorStats(stats, resid, NULL, mask, maskValue)) { + psError(PS_ERR_UNKNOWN, false, "Failure to compute statistics on the resid vector.\n"); + psFree(resid); + psFree(fit); + if (!inMask) psFree (mask); + return false; + } + + double meanValue = psStatsGetValue (stats, meanOption); + double stdevValue = psStatsGetValue (stats, stdevOption); + + psTrace("psLib.imageops", 5, "Mean is %f\n", meanValue); + psTrace("psLib.imageops", 5, "Stdev is %f\n", stdevValue); + psF32 minClipValue = -minClipSigma*stdevValue; + psF32 maxClipValue = +maxClipSigma*stdevValue; + + // set mask if pts are not valid + // we are masking out any point which is out of range + // recovery is not allowed with this scheme + for (psS32 i = 0; i < resid->n; i++) { + // XXX this prevents recovery of previously masked values + if (mask->data.PS_TYPE_VECTOR_MASK_DATA[i] & maskValue) { + continue; + } + + if ((resid->data.F32[i] - meanValue > maxClipValue) || (resid->data.F32[i] - meanValue < minClipValue)) { + psTrace("psLib.imageops", 6, "Masking element %d : %f vs %f : resid is %f\n", i, f->data.F32[i], fit->data.F32[i], resid->data.F32[i]); + mask->data.PS_TYPE_VECTOR_MASK_DATA[i] |= 0x01; + continue; + } + Nkeep++; + } + + // We should probably exit this loop if no new elements were masked since the fit won't + // change. + psTrace("psLib.imageops", 6, "keeping %d of %ld pts for fit\n", Nkeep, x->n); + stats->clippedNvalues = Nkeep; + psFree(fit); + } + + // Free local temporary variables + psFree(resid); + if (!inMask) psFree (mask); + *pGoodFit = true; // XXX probably don't need to set this (set by psImageMapFit) + return true; +} + +// CZW: 2014-10-09 +// Sparse versions of MapFit and MapFitClip that assume the matrices are not filled. +bool psImageMapFitSparse(bool *pGoodFit, psImageMap *map, const psVector *mask, psVectorMaskType maskValue, + const psVector *x, const psVector *y, const psVector *f, const psVector *df) +{ + // XXX Add Asserts + + *pGoodFit = false; + + // dimensions of the output map image + int Nx = map->binning->nXruff; + int Ny = map->binning->nYruff; + + + // no spatial information, just calculate mean & stdev + if ((Nx == 1) && (Ny == 1)) { + psStatsInit(map->stats); + + // the user has supplied one of various stats option pairs, + psStatsOptions mean = psStatsMeanOption(map->stats->options); + psStatsOptions stdev = psStatsStdevOption(map->stats->options); + if (!psStatsSingleOption(mean)) { + psError(PS_ERR_UNKNOWN, true, "no valid mean stats option selected"); + return false; + } + if (!psStatsSingleOption(stdev)) { + psError(PS_ERR_UNKNOWN, true, "no valid stdev stats option selected"); + return false; + } + + // XXX does ROBUST_MEDIAN work with weight? + if (!psVectorStats(map->stats, f, NULL, mask, maskValue)) { + psError(PS_ERR_UNKNOWN, false, "failure to measure stats"); + return false; + } + + map->map->data.F32[0][0] = psStatsGetValue(map->stats, mean); + map->error->data.F32[0][0] = psStatsGetValue(map->stats, stdev); + if (isfinite(map->map->data.F32[0][0]) && isfinite( map->error->data.F32[0][0])) { + *pGoodFit = true; + } + return true; + } + + if (Nx == 1) { + bool status; + status = psImageMapFit1DinY (pGoodFit, map, mask, maskValue, x, y, f, df); + return status; + } + if (Ny == 1) { + bool status; + status = psImageMapFit1DinX (pGoodFit, map, mask, maskValue, x, y, f, df); + return status; + } + + // set up the redirection table so we can use sA[-1][-1], etc + // XXX psKernel does this for you --- PAP. + float SAm[3][3], *SAv[3], **sA; + + for (int i = 0; i < 3; i++) { + SAv[i] = SAm[i] + 1; + } + sA = SAv + 1; + + // elements of the matrix equation Ax = B; we are solving for the vector x + // psImage *A = psImageAlloc (Nx*Ny, Nx*Ny, PS_TYPE_F32); + // psVector *B = psVectorAlloc (Nx*Ny, PS_TYPE_F32); + + // psImageInit (A, 0.0); + // psVectorInit (B, 0.0); + + // CZW: call to psSparseAlloc + // It should match old A, and each element of that should only touch four others. + psSparse *Asparse = psSparseAlloc(Nx * Ny, 4 * Nx * Ny); + + // we are looping over the Nx,Ny image map elements; + // the matrix equation contains Nx*Ny rows and columns + + for (int n = 0; n < Nx; n++) { + for (int m = 0; m < Ny; m++) { + // define & init summing variables + float rx_rx_ry_ry = 0; + float rx_rx_dy_ry = 0; + float dx_rx_ry_ry = 0; + float dx_rx_dy_ry = 0; + float fi_rx_ry = 0; + float rx_rx_py_py = 0; + float rx_rx_qy_py = 0; + float dx_rx_py_py = 0; + float dx_rx_qy_py = 0; + float fi_rx_py = 0; + float px_px_ry_ry = 0; + float px_px_dy_ry = 0; + float qx_px_ry_ry = 0; + float qx_px_dy_ry = 0; + float fi_px_ry = 0; + float px_px_py_py = 0; + float px_px_qy_py = 0; + float qx_px_py_py = 0; + float qx_px_qy_py = 0; + float fi_px_py = 0; + + // generate the sums for the fitting matrix element I,J + // I = n + nX*m + // J = (n + jn) + nX*(m + jm) + for (int i = 0; i < x->n; i++) { + + if (mask && (mask->data.PS_TYPE_VECTOR_MASK_DATA[i] & maskValue)) continue; + + // base coordinate offset for this point (x,y) relative to this map element (n,m) + float dx = psImageBinningGetRuffX (map->binning, x->data.F32[i]) - (n + 0.5); + float dy = psImageBinningGetRuffY (map->binning, y->data.F32[i]) - (m + 0.5); + + // edge cases to include: + bool edgeX = false; + edgeX |= ((n == 1) && (dx < -1.0)); + edgeX |= ((n == Nx - 2) && (dx > +1.0)); + + bool edgeY = false; + edgeY |= ((m == 1) && (dy < -1.0)); + edgeY |= ((m == Ny - 2) && (dy > +1.0)); + + // skip points outside of 2x2 grid centered on n,m: + if (!edgeX && (fabs(dx) > 1.0)) continue; + if (!edgeY && (fabs(dy) > 1.0)) continue; + + // related offset values + float rx = 1.0 - dx; + float ry = 1.0 - dy; + float px = 1.0 + dx; + float py = 1.0 + dy; + float qx = -dx; + float qy = -dy; + + // data value & weight for this point + float fi = f->data.F32[i]; + if (!isfinite(fi)) continue; + + float wt = 1.0; + if (df != NULL) { + if (df->data.F32[i] == 0.0) { + wt = 0.0; + } else { + if (!isfinite(df->data.F32[i])) continue; + wt = 1.0 / PS_SQR(df->data.F32[i]); // XXX test for dz == NULL or dz_i = 0 + } + } + + // sum the appropriate elements for the different quadrants + + int Qx = (dx >= 0) ? 1 : 0; + if (n == 0) Qx = 1; + if (n == Nx - 1) Qx = 0; + + int Qy = (dy >= 0) ? 1 : 0; + if (m == 0) Qy = 1; + if (m == Ny - 1) Qy = 0; + + assert (isfinite(fi)); + assert (isfinite(wt)); + assert (isfinite(rx)); + assert (isfinite(ry)); + + // points at offset 1,1 + if ((Qx == 1) && (Qy == 1)) { + rx_rx_ry_ry += rx*rx*ry*ry*wt; + rx_rx_dy_ry += rx*rx*dy*ry*wt; + dx_rx_ry_ry += dx*rx*ry*ry*wt; + dx_rx_dy_ry += dx*rx*dy*ry*wt; + fi_rx_ry += fi*rx*ry*wt; + } + // points at offset 1,0 + if ((Qx == 1) && (Qy == 0)) { + rx_rx_py_py += rx*rx*py*py*wt; + rx_rx_qy_py += rx*rx*qy*py*wt; + dx_rx_py_py += dx*rx*py*py*wt; + dx_rx_qy_py += dx*rx*qy*py*wt; + fi_rx_py += fi*rx*py*wt; + } + // points at offset 0,1 + if ((Qx == 0) && (Qy == 1)) { + px_px_ry_ry += px*px*ry*ry*wt; + px_px_dy_ry += px*px*dy*ry*wt; + qx_px_ry_ry += qx*px*ry*ry*wt; + qx_px_dy_ry += qx*px*dy*ry*wt; + fi_px_ry += fi*px*ry*wt; + } + // points at offset 0,0 + if ((Qx == 0) && (Qy == 0)) { + px_px_py_py += px*px*py*py*wt; + px_px_qy_py += px*px*qy*py*wt; + qx_px_py_py += qx*px*py*py*wt; + qx_px_qy_py += qx*px*qy*py*wt; + fi_px_py += fi*px*py*wt; + } + } + + // the chi-square derivatives have elements of the form g(n+jn,m+jm)*A(jn,jm), + // jn,jm = -1 to +1. Convert the sums above into the correct coefficients + sA[-1][-1] = qx_px_qy_py; + sA[-1][ 0] = qx_px_ry_ry + qx_px_py_py; + sA[-1][+1] = qx_px_dy_ry; + sA[ 0][-1] = rx_rx_qy_py + px_px_qy_py; + sA[ 0][ 0] = rx_rx_ry_ry + px_px_ry_ry + rx_rx_py_py + px_px_py_py; + sA[ 0][+1] = rx_rx_dy_ry + px_px_dy_ry; + sA[+1][-1] = dx_rx_qy_py; + sA[+1][ 0] = dx_rx_ry_ry + dx_rx_py_py; + sA[+1][+1] = dx_rx_dy_ry; + + // I[ 0][ 0] = index for this n,m element: + int I = n + Nx * m; + // B->data.F32[I] = fi_rx_ry + fi_rx_py + fi_px_ry + fi_px_py; + // CZW: call to psSparseVector Element + if (fi_rx_ry + fi_rx_py + fi_px_ry + fi_px_py == 0.0) { + psSparseVectorElement(Asparse, I, 1.0); + } + else { + psSparseVectorElement(Asparse, I, fi_rx_ry + fi_rx_py + fi_px_ry + fi_px_py); + } + + // printf("ADDING: %d %g \n",I, fi_rx_ry + fi_rx_py + fi_px_ry + fi_px_py); + // insert these values into their corresponding locations in A, B + for (int jn = -1; jn <= +1; jn++) { + if (n + jn < 0) continue; + if (n + jn >= Nx) continue; + for (int jm = -1; jm <= +1; jm++) { + if (m + jm < 0) continue; + if (m + jm >= Ny) continue; + int J = (n + jn) + Nx * (m + jm); + // printf("A: %d %d %g\n",J,I,sA[jn][jm]); + // A->data.F32[J][I] = sA[jn][jm]; + // CZW: call to psSparseMatrixElement + if (J < I) { continue; } + psSparseMatrixElement(Asparse,J,I,sA[jn][jm]); // Ensure J < I? + + } + } + } + } + + // test for empty diagonal elements (unconstained cells), mark, and set pivots to 1.0 + // CZW: I'm not totally sure how to check these in the sparse context. + // Iterate over all ii pairs, and manually check the structure? +#if (0) + psVector *Empty = psVectorAlloc (Nx*Ny, PS_TYPE_S8); + psVectorInit (Empty, 0); + for (int i = 0; i < Nx*Ny; i++) { + if (A->data.F32[i][i] == 0.0) { + Empty->data.S8[i] = 1; + for (int j = 0; j < Nx*Ny; j++) { + A->data.F32[i][j] = 0.0; + A->data.F32[j][i] = 0.0; + } + A->data.F32[i][i] = 1.0; + B->data.F32[i] = 0.0; + } + } +#endif + + // CZW: call to psSparseSolve + psVector *solution = psVectorAlloc(Nx*Ny, PS_TYPE_F32); + psSparseConstraint Constraint; + Constraint.paramDelta = 1e-3; + Constraint.paramMin = -1e5; + Constraint.paramMax = 1e5; + solution = psSparseSolve(solution, Constraint, Asparse, 1000); + if (!solution) { + psFree(solution); + psFree(Asparse); + return(false); + } + +#if (0) + // CZW: This is a continuation of the above information + // set bad values to NaN + for (int i = 0; i < Nx*Ny; i++) { + if (Empty->data.S8[i]) { + B->data.F32[i] = NAN; + A->data.F32[i][i] = 0; + } + } +#endif + + for (int n = 0; n < Nx; n++) { + for (int m = 0; m < Ny; m++) { + int I = n + Nx * m; + map->map->data.F32[m][n] = solution->data.F32[I]; + map->error->data.F32[m][n] = NAN; // sqrt(A->data.F32[I][I]); // CZW: fix this to be a real error. + } + } + + // psFree (A); + // psFree (B); + // psFree (Empty); + // CZW: free things + psFree(solution); + psFree(Asparse); + + *pGoodFit = true; + return true; +} + +// measure residuals on each pass and clip outliers based on stats +bool psImageMapClipFitSparse(bool *pGoodFit, psImageMap *map, psStats *stats, psVector *inMask, psVectorMaskType maskValue, + const psVector *x, const psVector *y, const psVector *f, const psVector *df) +{ + // XXX add in full PS_ASSERTS + psAssert(map, "impossible"); + psAssert(stats, "impossible"); + psAssert(x, "impossible"); + psAssert(y, "impossible"); + psAssert(f, "impossible"); + + *pGoodFit = false; + + // the user supplies one of various stats option pairs, + // determine the desired mean and stdev STATS options: + psStatsOptions meanOption = psStatsMeanOption(stats->options); + psStatsOptions stdevOption = psStatsStdevOption(stats->options); + if (!psStatsSingleOption(meanOption)) { + psError(PS_ERR_UNKNOWN, true, "no valid mean stats option selected"); + return false; + } + if (!psStatsSingleOption(stdevOption)) { + psError(PS_ERR_UNKNOWN, true, "no valid stdev stats option selected"); + return false; + } + + // clipping range defined by min and max and/or clipSigma + psF32 minClipSigma; + psF32 maxClipSigma; + if (isfinite(stats->max)) { + maxClipSigma = fabs(stats->max); + } else { + maxClipSigma = fabs(stats->clipSigma); + } + if (isfinite(stats->min)) { + minClipSigma = fabs(stats->min); + } else { + minClipSigma = fabs(stats->clipSigma); + } + + psVector *mask = inMask; + if (!inMask) { + mask = psVectorAlloc (x->n, PS_TYPE_VECTOR_MASK); + psVectorInit (mask, 0); + } + + // vector to store residuals + psVector *resid = psVectorAlloc(f->n, PS_TYPE_F32); + + psTrace("psLib.imageops", 4, "stats->clipIter is %d\n", stats->clipIter); + psTrace("psLib.imageops", 4, "(minClipSigma, maxClipSigma) is (%.2f, %.2f)\n", minClipSigma, maxClipSigma); + + for (psS32 N = 0; N < stats->clipIter; N++) { + psTrace("psLib.imageops", 6, "Loop iteration %d. Calling psImageMapFit()\n", N); + psS32 Nkeep = 0; + if (!psImageMapFitSparse(pGoodFit, map, mask, maskValue, x, y, f, df)) { psError(PS_ERR_UNKNOWN, false, "Could not fit image map.\n"); psFree(resid); Index: branches/eam_branches/ipp-20140904/psLib/src/imageops/psImageMapFit.h =================================================================== --- branches/eam_branches/ipp-20140904/psLib/src/imageops/psImageMapFit.h (revision 37803) +++ branches/eam_branches/ipp-20140904/psLib/src/imageops/psImageMapFit.h (revision 37804) @@ -20,4 +20,28 @@ // fit the image map to a set of points bool psImageMapClipFit(bool *pGoodFit, + psImageMap *map, + psStats *stats, + psVector *mask, // WARNING: Mask is modified! + psVectorMaskType maskValue, + const psVector *x, + const psVector *y, + const psVector *f, + const psVector *df + ); + + +// fit the image map to a set of points, using sparse matrix tools +bool psImageMapFitSparse(bool *pGoodFit, + psImageMap *map, + const psVector *mask, + psVectorMaskType maskValue, // + const psVector *x, + const psVector *y, + const psVector *f, + const psVector *df + ); + +// fit the image map to a set of points using sparse matrix tools +bool psImageMapClipFitSparse(bool *pGoodFit, psImageMap *map, psStats *stats, Index: branches/eam_branches/ipp-20140904/psLib/src/math/psMinimizePolyFit.c =================================================================== --- branches/eam_branches/ipp-20140904/psLib/src/math/psMinimizePolyFit.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psLib/src/math/psMinimizePolyFit.c (revision 37804) @@ -700,19 +700,27 @@ // Define values that may be used by PS_POLYNOMIAL_ORD form. + // these scaling values put the dynamic range of the data into some vaguely sensible location bool scale = false; + bool median_zero = false; double median = 0.0; double sigma = 1.0; - psVector *sorted = NULL; - psVector *z64 = NULL; - + switch (poly->type) { case PS_POLYNOMIAL_ORD: if ((f64->n < 10000)&&(poly->nX > 1)) { scale = true; - sorted = psVectorSort(NULL,x64); + + // generate a subset of the unmasked values: + psVector *tmp = psVectorAllocEmpty (x64->n, PS_TYPE_F64); + for (int itmp = 0; itmp < x64->n; itmp++) { + if (mask && (mask->data.PS_TYPE_VECTOR_MASK_DATA[itmp] && maskValue)) continue; + psVectorAppend (tmp, x64->data.F64[itmp]); + } + psVector *sorted = psVectorSort(NULL,tmp); median = sorted->data.F64[sorted->n / 2]; // CZW: I'm not bothering to scale this because it doesn't really matter. sigma = (sorted->data.F64[3 * sorted->n / 4] - sorted->data.F64[sorted->n / 4]); psFree(sorted); + psFree(tmp); if ((!isfinite(median))|| @@ -721,9 +729,16 @@ sigma = 1.0; } - if ((median == 0.0)&&(sigma == 1.0)) { - scale = false; + if (fabs(median) < 1e-10) { + // CZW 2014-10-21: This median is small and close to zero. This can cause issues with the + // scaling, + median_zero = true; + if ((sigma == 1.0)||(fabs(sigma) <= 1e-10)) { + // Don't bother scaling if sigma is unity (it's already scaled) or if the sigma calculation has gone wrong. + scale = false; + } } + // I can't see a way to not clobber x if it's already F64, so make a copy.x - z64 = psVectorCopy(NULL,x64,PS_TYPE_F64); + psVector *z64 = psVectorCopy(NULL,x64,PS_TYPE_F64); psBinaryOp(z64,z64,"-",psScalarAlloc(median,PS_TYPE_F64)); psBinaryOp(z64,z64,"/",psScalarAlloc(sigma,PS_TYPE_F64)); @@ -733,5 +748,4 @@ #endif - result = VectorFitPolynomial1DOrd(poly, mask, maskValue, f64, fErr64, z64); psFree(z64); // Done with this. @@ -758,22 +772,28 @@ #endif } + for (psS32 i = 0; i <= poly->nX; i++) { poly->coeff[i] = 0.0; poly->coeffErr[i] = 0.0; - - for (psS32 j = 0; j <= poly->nX; j++) { + if (median_zero) { // If the median is zero, the obtained solution just needs to be scaled by the sigma values. + poly->coeff[i] = Zcoeff[i] * pow(1.0 / sigma,i); + poly->coeffErr[i] = ZcoeffErr[i] * pow(1.0 / sigma,i); + } + else { // Otherwise, do the correct transformations by expanding the (x-m)/s terms. + for (psS32 j = 0; j <= poly->nX; j++) { #if (CZW) - printf(" %d %d %f %f %f %f => %f\n", - i,j,Zcoeff[j], - pow(1.0 / sigma,j) * pow(-1,j - i), - pow(median,j - i), - 1.0 * psBinomialCoeff(j,i), - Zcoeff[j] * pow(1.0 / sigma,j) * pow(-1,j -i) * pow(median,j - i) * 1.0 * psBinomialCoeff(j,i) - ); + printf(" %d %d %f %f %f %f => %f\n", + i,j,Zcoeff[j], + pow(1.0 / sigma,j) * pow(-1,j - i), + pow(median,j - i), + 1.0 * psBinomialCoeff(j,i), + Zcoeff[j] * pow(1.0 / sigma,j) * pow(-1,j -i) * pow(median,j - i) * 1.0 * psBinomialCoeff(j,i) + ); #endif - poly->coeff[i] += Zcoeff[j] * pow(1.0 / sigma,j) * pow(-1,j - i) * pow(median,j - i) * psBinomialCoeff(j,i); - poly->coeffErr[i] += pow(ZcoeffErr[j] * pow(1.0 / sigma,j) * pow(-1,j - i) * pow(median,j - 1) * psBinomialCoeff(j,i),2); + poly->coeff[i] += Zcoeff[j] * pow(1.0 / sigma,j) * pow(-1,j - i) * pow(median,j - i) * psBinomialCoeff(j,i); + poly->coeffErr[i] += pow(ZcoeffErr[j] * pow(1.0 / sigma,j) * pow(-1,j - i) * pow(median,j - 1) * psBinomialCoeff(j,i),2); + } + poly->coeffErr[i] = sqrt(poly->coeffErr[i]); } - poly->coeffErr[i] = sqrt(poly->coeffErr[i]); #if (CZW) printf("poly1d: unscaled parameters: %d %f %f\n", @@ -784,5 +804,5 @@ psFree(ZcoeffErr); - } + } // End scaling block. break; Index: branches/eam_branches/ipp-20140904/psLib/src/mathtypes/psImage.c =================================================================== --- branches/eam_branches/ipp-20140904/psLib/src/mathtypes/psImage.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psLib/src/mathtypes/psImage.c (revision 37804) @@ -91,5 +91,5 @@ } - long numBytes = numRows * numCols * elementSize; + size_t numBytes = (size_t) numRows * (size_t) numCols * (size_t) elementSize; psImage* image = (psImage* ) p_psAlloc(file, lineno, func, sizeof(psImage)); Index: branches/eam_branches/ipp-20140904/psModules/src/concepts/pmConcepts.c =================================================================== --- branches/eam_branches/ipp-20140904/psModules/src/concepts/pmConcepts.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psModules/src/concepts/pmConcepts.c (revision 37804) @@ -312,4 +312,7 @@ conceptRegisterEnum("FPA.TIMESYS", "Time system", p_pmConceptParse_TIMESYS, p_pmConceptFormat_TIMESYS, p_pmConceptCopy_TIMESYS, false, PM_FPA_LEVEL_FPA); conceptRegisterTime("FPA.TIME", "Time of exposure", false, PM_FPA_LEVEL_FPA); + + conceptRegisterStr("FPA.SHUTOUTC", "Time of exposure", NULL,NULL,NULL,false, PM_FPA_LEVEL_FPA); + conceptRegisterF32("FPA.TEMP", "Temperature of focal plane", NULL, NULL, NULL, false, PM_FPA_LEVEL_FPA); conceptRegisterF32("FPA.M1X", "Primary Mirror X Position", NULL, NULL, NULL, false, PM_FPA_LEVEL_FPA); Index: branches/eam_branches/ipp-20140904/psModules/src/imcombine/pmSubtraction.c =================================================================== --- branches/eam_branches/ipp-20140904/psModules/src/imcombine/pmSubtraction.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psModules/src/imcombine/pmSubtraction.c (revision 37804) @@ -995,4 +995,19 @@ } + + if (0) { + FILE *f = fopen ("vector.dat", "w"); + psAssert (f, "should not fail to open test file"); + + for (int qq = 0; qq < match->chisq->n; qq++) { + fprintf (f, "%f %f : %f %f : %d\n", + match->chisq->data.F32[qq], match->fluxes->data.F32[qq], + logchisq->data.F32[qq], logflux->data.F32[qq], + match->stampMask->data.PS_TYPE_VECTOR_MASK_DATA[qq]); + } + fclose (f); + } + + bool result = psVectorClipFitPolynomial1D(model, stats, match->stampMask, 0xff, logchisq, NULL, logflux); if (!result) { Index: branches/eam_branches/ipp-20140904/psModules/src/objects/pmSourceExtendedPars.c =================================================================== --- branches/eam_branches/ipp-20140904/psModules/src/objects/pmSourceExtendedPars.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psModules/src/objects/pmSourceExtendedPars.c (revision 37804) @@ -226,4 +226,10 @@ pars->petrosianR50 = NAN; pars->petrosianR50Err = NAN; + pars->ghalfLightRadius = NAN; + pars->gRT = NAN; + pars->gRA = NAN; + pars->gS2 = NAN; + pars->gA = NAN; + pars->gbumpy = NAN; return pars; } Index: branches/eam_branches/ipp-20140904/psModules/src/objects/pmSourceExtendedPars.h =================================================================== --- branches/eam_branches/ipp-20140904/psModules/src/objects/pmSourceExtendedPars.h (revision 37803) +++ branches/eam_branches/ipp-20140904/psModules/src/objects/pmSourceExtendedPars.h (revision 37804) @@ -65,4 +65,10 @@ float petrosianR50Err; float petrosianFill; + float ghalfLightRadius; + float gRT; // total residual (from eliptically symmetric model) Simard 2002 + float gRA; // assymetric residual (from eliptically symmetric model) Simard 2002 + float gS2; // (un) smoothness (Simard, 2011 Cheng 2011) + float gA; // assymetry index (Gyory & Bell 2010) + float gbumpy; // Blakeslee bumpiness } pmSourceExtendedPars; Index: branches/eam_branches/ipp-20140904/psModules/src/objects/pmSourceIO_CMF.c.in =================================================================== --- branches/eam_branches/ipp-20140904/psModules/src/objects/pmSourceIO_CMF.c.in (revision 37803) +++ branches/eam_branches/ipp-20140904/psModules/src/objects/pmSourceIO_CMF.c.in (revision 37804) @@ -319,5 +319,6 @@ // define PSF model type - int modelType = pmModelClassGetType ("PS_MODEL_GAUSS"); + int defaultModelType = pmModelClassGetType ("PS_MODEL_GAUSS"); + int modelType = -1; // if header does not define the model, default to a gaussian @@ -326,5 +327,9 @@ modelType = pmModelClassGetType (PSF_NAME); } - assert (modelType > -1); + // work around bug in psphotFullForce + if (modelType < 0) { + modelType = defaultModelType; + } + // assert (modelType > -1); // do we expect to find lensing parameters? @@ -456,4 +461,8 @@ @>PS1_V4@ source->lensingOBJ->shear->e2 = psMetadataLookupF32 (&status, row, "LENS_E2_SH_OBJ"); } + + @>PS1_V4@ source->chipNum = psMetadataLookupS16 (&status, row, "SRC_CHIP_NUM"); + @>PS1_V4@ source->chipX = psMetadataLookupS16 (&status, row, "SRC_CHIP_X"); + @>PS1_V4@ source->chipY = psMetadataLookupS16 (&status, row, "SRC_CHIP_Y"); if (haveLensPSF) { @@ -639,4 +648,10 @@ psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_90_ERR", PS_DATA_F32, "Petrosian R90 Error (pix)", extpars->petrosianR90Err); psMetadataAdd (row, PS_LIST_TAIL, "PETRO_FILL", PS_DATA_F32, "Petrosian Fill Factor", extpars->petrosianFill); + @>PS1_SV3@ psMetadataAddF32 (row, PS_LIST_TAIL, "HALF_LIGHT_RADIUS", 0, "Half Light Radius", extpars->ghalfLightRadius); + @>PS1_SV3@ psMetadataAddF32 (row, PS_LIST_TAIL, "G_RT", 0, "total residual Sersic Model", extpars->gRT); + @>PS1_SV3@ psMetadataAddF32 (row, PS_LIST_TAIL, "G_RA", 0, "assymetric residual Sersic Model", extpars->gRA); + @>PS1_SV3@ psMetadataAddF32 (row, PS_LIST_TAIL, "G_S2", 0, "(un)smoothness paramter residual Sersic Model", extpars->gS2); + @>PS1_SV3@ psMetadataAddF32 (row, PS_LIST_TAIL, "G_A", 0, "assymetry index", extpars->gA); + @>PS1_SV3@ psMetadataAddF32 (row, PS_LIST_TAIL, "G_BUMPY", 0, "bumpiness", extpars->gbumpy); } else { psMetadataAdd (row, PS_LIST_TAIL, "PETRO_MAG", PS_DATA_F32, "Petrosian Magnitude", NAN); @@ -649,4 +664,10 @@ psMetadataAdd (row, PS_LIST_TAIL, "PETRO_RADIUS_90_ERR", PS_DATA_F32, "Petrosian R90 Error (pix)",NAN); psMetadataAdd (row, PS_LIST_TAIL, "PETRO_FILL", PS_DATA_F32, "Petrosian Fill Factor", NAN); + @>PS1_SV3@ psMetadataAddF32 (row, PS_LIST_TAIL, "HALF_LIGHT_RADIUS", 0, "Half Light Radius", NAN); + @>PS1_SV3@ psMetadataAddF32 (row, PS_LIST_TAIL, "G_RT", 0, "total residual Sersic Model", NAN); + @>PS1_SV3@ psMetadataAddF32 (row, PS_LIST_TAIL, "G_RA", 0, "assymetric residual Sersic Model", NAN); + @>PS1_SV3@ psMetadataAddF32 (row, PS_LIST_TAIL, "G_S2", 0, "(un)smoothness paramter residual Sersic Model", NAN); + @>PS1_SV3@ psMetadataAddF32 (row, PS_LIST_TAIL, "G_A", 0, "assymetry index", NAN); + @>PS1_SV3@ psMetadataAddF32 (row, PS_LIST_TAIL, "G_BUMPY", 0, "bumpiness", NAN); } } @@ -848,4 +869,10 @@ extpars->petrosianR90Err = psMetadataLookupF32(&status, row, "PETRO_RADIUS_90_ERR"); extpars->petrosianFill = psMetadataLookupF32(&status, row, "PETRO_FILL"); + extpars->ghalfLightRadius = psMetadataLookupF32(&status, row, "HALF_LIGHT_RADIUS"); + extpars->gRT = psMetadataLookupF32(&status, row, "G_RT"); + extpars->gRA = psMetadataLookupF32(&status, row, "G_RA"); + extpars->gS2 = psMetadataLookupF32(&status, row, "G_S2"); + extpars->gA = psMetadataLookupF32(&status, row, "G_A"); + extpars->gbumpy = psMetadataLookupF32(&status, row, "G_BUMPY"); psVector *radSB = psMetadataLookupVector(&status, row, "PROF_SB"); @@ -995,4 +1022,5 @@ psMetadataAddF32 (row, PS_LIST_TAIL, "X_EXT_SIG", 0, "Sigma in EXT x coordinate", xErr); psMetadataAddF32 (row, PS_LIST_TAIL, "Y_EXT_SIG", 0, "Sigma in EXT y coordinate", yErr); + @>PS1_SV3@ psMetadataAddF32 (row, PS_LIST_TAIL, "SKY_EXT", 0, "EXT fit sky value", PAR ? PAR[PM_PAR_SKY] : NAN); @>PS1_DV2,>PS1_SV3@ psMetadataAddF32 (row, PS_LIST_TAIL, "RA_EXT", 0, "EXT model ra coordinate", raPos); @>PS1_DV2,>PS1_SV3@ psMetadataAddF32 (row, PS_LIST_TAIL, "DEC_EXT", 0, "EXT model dec coordinate", decPos); Index: branches/eam_branches/ipp-20140904/psastro/src/psastro.h =================================================================== --- branches/eam_branches/ipp-20140904/psastro/src/psastro.h (revision 37803) +++ branches/eam_branches/ipp-20140904/psastro/src/psastro.h (revision 37804) @@ -154,7 +154,8 @@ bool psastroZeroPoint (pmConfig *config); -psVector *psastroZeroPointReadoutAccum(psVector *dMag, pmReadout *readout, float exptime); +psVector *psastroZeroPointReadoutAccum(psVector *dMag, pmReadout *readout, float exptime, float refminMag, float refmaxMag); bool psastroZeroPointAnalysis (psMetadata *header, psVector *dMag, float zeropt, psMetadata *recipe); bool psastroZeroPointFromRecipe (float *zeropt, float *exptime, float *ghostMaxMag, pmFPA *fpa, psMetadata *recipe); +bool psastroZeroPointRefMagLimitFromRecipe (float *refminMag, float *refmaxMag, pmFPA *fpa, psMetadata *recipe); psStats *psastroStatsPercentile (psVector *myVector, psMetadata *recipe); Index: branches/eam_branches/ipp-20140904/psastro/src/psastroChooseRefstars.c =================================================================== --- branches/eam_branches/ipp-20140904/psastro/src/psastroChooseRefstars.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psastro/src/psastroChooseRefstars.c (revision 37804) @@ -66,4 +66,9 @@ psF32 clampMagMin = psMetadataLookupF32 (&status, recipe, "REFSTAR_CLAMP_MAG_MIN"); if (!status) clampMagMin = -5.0; + //MEH adding max option as well -- really should not be hardcoded.. + psF32 clampMagMax = psMetadataLookupF32 (&status, recipe, "REFSTAR_CLAMP_MAG_MAX"); + if (!status) clampMagMax = 30.0; + psWarning("Will skip/clamp refstar with brighter/fainter magnitude %f/%f\n",clampMagMin,clampMagMax); + // de-activate all files except PSASTRO.REFSTARS @@ -109,6 +114,16 @@ if (ref->Mag < clampMagMin) { - psWarning("Skipping refstar with abusrd magnitude %f ra: %f dec: %f\n", - ref->Mag, RAD_TO_DEG(ref->sky->r), RAD_TO_DEG(ref->sky->d)); + //psWarning("Skipping refstar with abusrd magnitude %f ra: %f dec: %f\n", + // ref->Mag, RAD_TO_DEG(ref->sky->r), RAD_TO_DEG(ref->sky->d)); + psTrace ("psastro", 6, "Skipping by clamp refstar with bright magnitude %f ra: %f dec: %f\n", + ref->Mag, RAD_TO_DEG(ref->sky->r), RAD_TO_DEG(ref->sky->d)); + goto skip; + } + //MEH option for max mag clamp also -- need to change psWarning (to log) to trace and above for bright.. + if (ref->Mag > clampMagMax) { + //psWarning("Skipping by clamp refstar with faint magnitude %f ra: %f dec: %f\n", + // ref->Mag, RAD_TO_DEG(ref->sky->r), RAD_TO_DEG(ref->sky->d)); + psTrace ("psastro", 6, "Skipping by clamp refstar with faint magnitude %f ra: %f dec: %f\n", + ref->Mag, RAD_TO_DEG(ref->sky->r), RAD_TO_DEG(ref->sky->d)); goto skip; } Index: branches/eam_branches/ipp-20140904/psastro/src/psastroZeroPoint.c =================================================================== --- branches/eam_branches/ipp-20140904/psastro/src/psastroZeroPoint.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psastro/src/psastroZeroPoint.c (revision 37804) @@ -20,4 +20,5 @@ pmCell *cell = NULL; pmReadout *readout = NULL; + float refminMag,refmaxMag; // select the current recipe @@ -59,6 +60,15 @@ } + // MEH -- RefMagLimit for ZP with photcode option (merr probably better but not an option) + // -- defualt is to use general config value (no hardcode..) + if (!psastroZeroPointRefMagLimitFromRecipe (&refminMag, &refmaxMag, fpa, recipe)) { + psLogMsg ("psastro", PS_LOG_INFO, "failed to load ref min,maxMag data from recipe"); + return false; + } + // if we measure the zero point by exposure, accumulate the dMag values here: psVector *dMag = NULL; + // MEH -- also the number of refcat sources used + int zpt_nref = 0; float fpaZP = 0.0; // Average zero point @@ -81,10 +91,16 @@ // calculate dMag for the matched stars - dMag = psastroZeroPointReadoutAccum (dMag, readout, exptime); - + //dMag = psastroZeroPointReadoutAccum (dMag, readout, exptime); + // MEH -- additional limits should be available for ZP if just robust mean and w/ stack skycal + dMag = psastroZeroPointReadoutAccum (dMag, readout, exptime, refminMag, refmaxMag); if (!byExposure) { // calculate dMag for the matched stars just for this readout (well, chip) psMetadata *header = psMetadataLookupMetadata (&status, readout->analysis, "PSASTRO.HEADER"); psastroZeroPointAnalysis (header, dMag, zeropt, recipe); + // MEH -- per readout N refcat, min/max mag -- added here, not part of median/edge process + zpt_nref += dMag->n; + psMetadataAddS32 (header, PS_LIST_TAIL, "ZPT_NREF", PS_META_REPLACE, "N refcat sources for zero point", dMag->n); + psMetadataAddF32 (header, PS_LIST_TAIL, "ZPT_MIN", PS_META_REPLACE, "min refcat mag for zero point", refminMag); + psMetadataAddF32 (header, PS_LIST_TAIL, "ZPT_MAX", PS_META_REPLACE, "max refcat mag for zero point", refmaxMag); psFree (dMag); dMag = NULL; @@ -95,4 +111,5 @@ numZP++; } + @@ -110,4 +127,6 @@ } psastroZeroPointAnalysis (header, dMag, zeropt, recipe); + //MEH -- all N refcat sources here + zpt_nref = dMag->n; psFree (dMag); dMag = NULL; @@ -141,9 +160,12 @@ // calculate dMag for the matched stars just for this readout (well, chip) psMetadata *header = psMetadataLookupMetadata (&status, readout->analysis, "PSASTRO.HEADER"); - psMetadataAddF32 (header, PS_LIST_TAIL, "ZPT_OBS", PS_META_REPLACE, "measured zero point", zptObs); psMetadataAddF32 (header, PS_LIST_TAIL, "ZPT_REF", PS_META_REPLACE, "reference zero point", zptRef); psMetadataAddF32 (header, PS_LIST_TAIL, "ZPT_ERR", PS_META_REPLACE, "error on zero point", zptErr); psMetadataAddF32 (header, PS_LIST_TAIL, "ZPT_OFF", PS_META_REPLACE, "zero point offset", zptOff); + //MEH -- per readout N refcat, MIN/MAX for consistency + psMetadataAddS32 (header, PS_LIST_TAIL, "ZPT_NREF", PS_META_REPLACE, "N refcat sources for zero point", zpt_nref); + psMetadataAddF32 (header, PS_LIST_TAIL, "ZPT_MIN", PS_META_REPLACE, "min refcat mag for zero point", refminMag); + psMetadataAddF32 (header, PS_LIST_TAIL, "ZPT_MAX", PS_META_REPLACE, "max refcat mag for zero point", refmaxMag); } } @@ -151,4 +173,14 @@ } + // MEH add N refcat, min/max mag to primary header + psMetadata *header = psMetadataLookupMetadata (&status, fpa->analysis, "PSASTRO.HEADER"); + if (!header) { + header = psMetadataAlloc (); + psMetadataAddMetadata (fpa->analysis, PS_LIST_TAIL, "PSASTRO.HEADER", PS_META_REPLACE, "psastro header stats", header); + psFree (header); + } + psMetadataAddS32 (header, PS_LIST_TAIL, "ZPT_NREF", PS_META_REPLACE, "Total N refcat sources for zero point", zpt_nref); + psMetadataAddF32 (header, PS_LIST_TAIL, "ZPT_MIN", PS_META_REPLACE, "min refcat mag for zero point", refminMag); + psMetadataAddF32 (header, PS_LIST_TAIL, "ZPT_MAX", PS_META_REPLACE, "max refcat mag for zero point", refmaxMag); psMetadataItem *item = psMetadataLookup(fpa->concepts, "FPA.ZP"); @@ -162,5 +194,5 @@ * accumulate the dMag values from this readout */ -psVector *psastroZeroPointReadoutAccum(psVector *dMag, pmReadout *readout, float exptime) { +psVector *psastroZeroPointReadoutAccum(psVector *dMag, pmReadout *readout, float exptime, float refminMag, float refmaxMag) { bool status; @@ -191,6 +223,9 @@ // XXX I should be applying the airmass term here or in psastroLoadRefstars // (ie, ref->Mag would be a prediction for a given star. - float value = ref->Mag - raw->Mag - 2.5*log10(exptime); - psVectorAppend (dMag, value); + // MEH -- only append if in allowed ref mag range -- may want merr but not available.. + if ( ref->Mag >= refminMag && ref->Mag <= refmaxMag ){ + float value = ref->Mag - raw->Mag - 2.5*log10(exptime); + psVectorAppend (dMag, value); + } } return dMag; @@ -490,5 +525,9 @@ continue; } - } + //MEH null is a pain.. so only log if set + psLogMsg ("psastro", PS_LOG_INFO, "found GHOST_MAX_MAG %f",*ghostMaxMag); + } + //MEH what zpt is set to + psLogMsg ("psastro", PS_LOG_INFO, "found ZEROPT %f",*zeropt); psFree (iter); return true; @@ -497,2 +536,75 @@ return false; } + + + +bool psastroZeroPointRefMagLimitFromRecipe (float *refminMag, float *refmaxMag, pmFPA *fpa, psMetadata *recipe) { + + bool status; + float tmpmin,tmpmax; + + // MEH -- duplicate psastroZeroPointFromRecipe for ZeroPoint adding RefMag limits with photcode option for extra restrictions + // must have non-photcode dependent setting -- use as default rather than a hardcoded value + *refminMag = psMetadataLookupF32 (&status, recipe, "REFSTAR.ZP.MIN.MAG"); + //if (!status) *refminMag = -5.0; + if (!status) { + psLogMsg ("psastro", PS_LOG_INFO, "RefMagLimit is missing REFSTAR.ZP.MIN.MAG"); + return false; + } + *refmaxMag = psMetadataLookupF32 (&status, recipe, "REFSTAR.ZP.MAX.MAG"); + //if (!status) *refmaxMag = 30.0; + if (!status) { + psLogMsg ("psastro", PS_LOG_INFO, "RefMagLimit is missing REFSTAR.ZP.MAX.MAG"); + return false; + } + + // select the filter; default to any photcode and mag limit otherwise + char *filter = psMetadataLookupStr (&status, fpa->concepts, "FPA.FILTERID"); + if (!status) ESCAPE ("RefMagLimit missing FPA.FILTER in concepts"); + + // need to select the PHOTCODE.DATA folder that matches our filter + psMetadataItem *item = psMetadataLookup (recipe, "PHOTCODE.DATA"); + if (!item) ESCAPE ("RefMagLimit PHOTCODE.DATA folders missing"); + if (item->type != PS_DATA_METADATA_MULTI) ESCAPE ("RefMagLimit PHOTCODE.DATA not a multi"); + + // PHOTCODE.DATA is a multi of metadata items + psListIterator *iter = psListIteratorAlloc(item->data.list, PS_LIST_HEAD, false); + + psMetadataItem *refItem = NULL; + while ((refItem = psListGetAndIncrement (iter))) { + if (refItem->type != PS_DATA_METADATA) ESCAPE ("RefMagLimit PHOTCODE.DATA entry is not a metadata folder"); + + char *refFilter = psMetadataLookupStr (&status, refItem->data.md, "FILTER"); + if (!status) { + // always seems to happen once, commented out in ZeroPoint above as well + //psLogMsg ("psastro", PS_LOG_INFO, "RefMagLimit PHOTCODE.DATA recipe folder is missing FILTER"); + continue; + } + + // does this entry match the current filter? + if (strcmp (refFilter, filter)) continue; + + psLogMsg ("psastro", PS_LOG_DETAIL, "RefMagLimit PHOTCODE.DATA found for filter %s", filter); + + // -- may set one or other or both -- no reason to log if found or often not + tmpmin = psMetadataLookupF32 (&status, refItem->data.md, "REFSTAR.ZP.MIN.MAG"); + if (status) { + *refminMag = tmpmin; + } + tmpmax = psMetadataLookupF32 (&status, refItem->data.md, "REFSTAR.ZP.MAX.MAG"); + if (status) { + *refmaxMag = tmpmax; + } + // once finds filter, should free and return true but multiple groups cause issue so just go through all groups now + //psFree (iter); + //return true; + } + + // log out set values + psLogMsg ("psastro", PS_LOG_INFO, "using REFSTAR.ZP.MIN.MAG %f",*refminMag); + psLogMsg ("psastro", PS_LOG_INFO, "using REFSTAR.ZP.MAX.MAG %f",*refmaxMag); + + psFree (iter); + return true; +} + Index: branches/eam_branches/ipp-20140904/psphot/src/Makefile.am =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/Makefile.am (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/Makefile.am (revision 37804) @@ -260,5 +260,6 @@ psphotSetNFrames.c \ psphotSourceMemory.c \ - psphotChipParams.c + psphotChipParams.c \ + psphotGalaxyParams.c # not currently used Index: branches/eam_branches/ipp-20140904/psphot/src/psphot.h =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphot.h (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphot.h (revision 37804) @@ -138,9 +138,13 @@ bool psphotExtendedSourceFits_Threaded (psThreadJob *job); +bool psphotGalaxyParams (pmConfig *config, const pmFPAview *view, const char *filerule); +bool psphotGalaxyParamsReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe); +bool psphotGalaxyParams_Threaded (psThreadJob *job); + bool psphotApResid (pmConfig *config, const pmFPAview *view, const char *filerule); bool psphotApResidReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe); bool psphotMagnitudes (pmConfig *config, const pmFPAview *view, const char *filerule); -bool psphotMagnitudesReadout(pmConfig *config, psMetadata *recipe, const pmFPAview *view, pmReadout *readout, psArray *sources, pmPSF *psf); +bool psphotMagnitudesReadout(pmConfig *config, psMetadata *recipe, const pmFPAview *view, pmReadout *readout, psArray *sources, pmPSF *psf, int index); bool psphotMagnitudes_Threaded (psThreadJob *job); @@ -168,4 +172,5 @@ // in psphotChoosePSF.c: bool psphotPSFstats (pmReadout *readout, pmPSF *psf); +bool psphotPSFstatsSources (pmReadout *readout, psArray *sources, pmPSF *psf); bool psphotMomentsStats (pmReadout *readout, psArray *sources); @@ -232,6 +237,6 @@ int psphotSaveImage (psMetadata *header, psImage *image, char *filename); bool psphotDumpConfig (pmConfig *config); -pmReadout *psphotSelectBackground (pmConfig *config, const pmFPAview *view); -pmReadout *psphotSelectBackgroundStdev (pmConfig *config, const pmFPAview *view); +pmReadout *psphotSelectBackground (pmConfig *config, const pmFPAview *view, int index); +pmReadout *psphotSelectBackgroundStdev (pmConfig *config, const pmFPAview *view, int index); // PSF / DBL / EXT evaluation functions Index: branches/eam_branches/ipp-20140904/psphot/src/psphotChooseAnalysisOptions.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotChooseAnalysisOptions.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotChooseAnalysisOptions.c (revision 37804) @@ -606,4 +606,10 @@ psLogMsg ("psphot.options", PS_LOG_WARN, "choose analysis options for %ld objects: %f sec\n", objects->n, psTimerMark ("psphot.options")); + psFree(exptime); + psFree(zeropt); + psFree(chips); + psFree(inputFilters); + psFree(petroFluxLim); + psFree(extFitFluxLim); psFree(toGal); Index: branches/eam_branches/ipp-20140904/psphot/src/psphotChoosePSF.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotChoosePSF.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotChoosePSF.c (revision 37804) @@ -1,3 +1,4 @@ # include "psphotInternal.h" +bool psphotPSFstatsSources (pmReadout *readout, psArray *sources, pmPSF *psf); // generate a PSF model for inputs without PSF models already loaded @@ -388,5 +389,7 @@ psFree (models); - if (!psphotPSFstats (readout, psf)) { + // if (!psphotPSFstats (readout, psf)) { + + if (!psphotPSFstatsSources (readout, sources, psf)) { psError(PSPHOT_ERR_PSF, false, "cannot measure PSF shape terms"); psFree(options); @@ -466,13 +469,13 @@ if (modelPSF->params->n > 7) { - psVectorAppend (psfExtra1, modelPSF->params->data.F32[7]); + psVectorAppend (psfExtra1, modelPSF->params->data.F32[7]); } if (modelPSF->params->n > 8) { - psVectorAppend (psfExtra2, modelPSF->params->data.F32[8]); + psVectorAppend (psfExtra2, modelPSF->params->data.F32[8]); } psFree (modelPSF); - // float fwhmtest = pmPSFtoFWHM(psf, xc, yc); - // fprintf (stderr, "fwhm: %f, %f : %f\n", FWHM_MAJOR, FWHM_MINOR, fwhmtest); + float fwhmtest = pmPSFtoFWHM(psf, xc, yc); + fprintf (stderr, "fwhm: %f, %f : %f\n", FWHM_MAJOR, FWHM_MINOR, fwhmtest); } } @@ -522,19 +525,172 @@ if (psfExtra1->n) { - if (!psVectorStats (stats, psfExtra1, NULL, NULL, 0)) { - psError(PS_ERR_UNKNOWN, false, "failure to measure stats for PSF EXTRA 1"); + if (!psVectorStats (stats, psfExtra1, NULL, NULL, 0)) { + psError(PS_ERR_UNKNOWN, false, "failure to measure stats for PSF EXTRA 1"); + goto escape; + } + psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "PSF_EXT1", PS_META_REPLACE, "PSF extra param 1", stats->sampleMean); + psLogMsg ("psphot", PS_LOG_DETAIL, "PSF extra parameter 1: %f +/- %f", stats->sampleMean, stats->sampleStdev); + } + + if (psfExtra2->n) { + if (!psVectorStats (stats, psfExtra2, NULL, NULL, 0)) { + psError(PS_ERR_UNKNOWN, false, "failure to measure stats for PSF EXTRA 2"); + goto escape; + } + psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "PSF_EXT2", PS_META_REPLACE, "PSF extra param 2", stats->sampleMean); + psLogMsg ("psphot", PS_LOG_DETAIL, "PSF extra parameter 2: %f +/- %f", stats->sampleMean, stats->sampleStdev); + } + + psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "ANGLE", PS_META_REPLACE, "PSF angle", axes.theta); + psMetadataAddS32 (readout->analysis, PS_LIST_TAIL, "NPSFSTAR", PS_META_REPLACE, "Number of stars used to make PSF", psf->nPSFstars); + + // psLogMsg ("psphot", PS_LOG_DETAIL, "PSF angle: %f, nstars: %d", axes.theta, psf->nPSFstars); + + char *psfModelName = pmModelClassGetName(psf->type); + psMetadataAddStr(readout->analysis, PS_LIST_TAIL, "PSFMODEL", PS_META_REPLACE, "PSF Model Name", psfModelName); + psMetadataAddBool(readout->analysis, PS_LIST_TAIL, "PSF_OK", PS_META_REPLACE, "Valid PSF Model?", true); + + int nParams = pmModelClassParameterCount(psf->type); + psMetadataAddS32(readout->analysis, PS_LIST_TAIL, "PSF_NPAR", PS_META_REPLACE, "Number of PSF parameters", nParams); + + psFree (fwhmMajor); + psFree (fwhmMinor); + psFree (psfExtra1); + psFree (psfExtra2); + psFree (stats); + return true; + +escape: + psFree (fwhmMajor); + psFree (fwhmMinor); + psFree (psfExtra1); + psFree (psfExtra2); + psFree (stats); + return false; +} + +// measure average parameters of the PSF model +bool psphotPSFstatsSources (pmReadout *readout, psArray *sources, pmPSF *psf) { + + psEllipseShape shape; + psEllipseAxes axes; + + PS_ASSERT_PTR_NON_NULL(readout, false); + PS_ASSERT_PTR_NON_NULL(psf, false); + + psImage *image = readout->image; + PS_ASSERT_PTR_NON_NULL(image, false); + + // XXX a minor hack: measure the values for a grid of points across the image, determine mean, UQ, LQ, etc: + psVector *fwhmMajor = psVectorAllocEmpty (100, PS_DATA_F32); + psVector *fwhmMinor = psVectorAllocEmpty (100, PS_DATA_F32); + psVector *psfExtra1 = psVectorAllocEmpty (100, PS_DATA_F32); + psVector *psfExtra2 = psVectorAllocEmpty (100, PS_DATA_F32); + + for (int i = 0; i < sources->n; i++) { + + pmSource *source = sources->data[i]; + + if (!(source->mode & PM_SOURCE_MODE_PSFSTAR)) continue; + + float xc = source->peak->xf; + float yc = source->peak->yf; + + // create modelPSF from this model + pmModel *modelPSF = pmModelFromPSFforXY (psf, xc, yc, 1.0); + if (!modelPSF) { + fprintf (stderr, "?"); + continue; + } + + // get the model full-width at half-max + float FWHM_MAJOR = 2*modelPSF->class->modelRadius (modelPSF->params, 0.5); + + // XXX make sure this is consistent with the re-definition of PM_PAR_SXX + shape.sx = modelPSF->params->data.F32[PM_PAR_SXX]; + shape.sy = modelPSF->params->data.F32[PM_PAR_SYY]; + shape.sxy = modelPSF->params->data.F32[PM_PAR_SXY]; + axes = psEllipseShapeToAxes (shape, 20.0); + + float FWHM_MINOR = FWHM_MAJOR * (axes.minor / axes.major); + if (!isfinite(FWHM_MAJOR) || !isfinite(FWHM_MINOR)) { + fprintf (stderr, "!"); + psFree (modelPSF); + continue; + } + psVectorAppend (fwhmMajor, FWHM_MAJOR); + psVectorAppend (fwhmMinor, FWHM_MINOR); + + if (modelPSF->params->n > 7) { + psVectorAppend (psfExtra1, modelPSF->params->data.F32[7]); + } + if (modelPSF->params->n > 8) { + psVectorAppend (psfExtra2, modelPSF->params->data.F32[8]); + } + psFree (modelPSF); + + float fwhmtest = pmPSFtoFWHM(psf, xc, yc); + fprintf (stderr, "fwhm: %f, %f : %f\n", FWHM_MAJOR, FWHM_MINOR, fwhmtest); + } + + psStats *stats = psStatsAlloc (PS_STAT_SAMPLE_MEAN | PS_STAT_SAMPLE_STDEV | PS_STAT_SAMPLE_QUARTILE); + if (!psVectorStats (stats, fwhmMajor, NULL, NULL, 0)) { + psError(PS_ERR_UNKNOWN, false, "failure to measure stats for FWHM MAJOR"); goto escape; - } - psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "PSF_EXT1", PS_META_REPLACE, "PSF extra param 1", stats->sampleMean); - psLogMsg ("psphot", PS_LOG_DETAIL, "PSF extra parameter 1: %f +/- %f", stats->sampleMean, stats->sampleStdev); + } + + psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FWHM_MAJ", PS_META_REPLACE, "PSF FWHM Major axis (mean)", stats->sampleMean); + psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MJ_SG", PS_META_REPLACE, "PSF FWHM Major axis (sigma)", stats->sampleStdev); + psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MJ_LQ", PS_META_REPLACE, "PSF FWHM Major axis (lower quartile)", stats->sampleLQ); + psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MJ_UQ", PS_META_REPLACE, "PSF FWHM Major axis (upper quartile)", stats->sampleUQ); + + float fwhmMaj = stats->sampleMean; // FWHM on major axis + + if (!psVectorStats (stats, fwhmMinor, NULL, NULL, 0)) { + psError(PS_ERR_UNKNOWN, false, "failure to measure stats for FWHM MINOR"); + goto escape; + } + psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FWHM_MIN", PS_META_REPLACE, "PSF FWHM Minor axis (mean)", stats->sampleMean); + psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MN_SG", PS_META_REPLACE, "PSF FWHM Minor axis (sigma)", stats->sampleStdev); + psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MN_LQ", PS_META_REPLACE, "PSF FWHM Minor axis (lower quartile)", stats->sampleLQ); + psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "FW_MN_UQ", PS_META_REPLACE, "PSF FWHM Minor axis (upper quartile)", stats->sampleUQ); + + float fwhmMin = stats->sampleMean; // FWHM on minor axis + + if (readout->parent) { + + // we now have 2 possible measurements of the seeing : the PSF model based version and + // the Moments based value we need to define a definitive "CHIP.SEEING" value. This + // value is used by the PSF-matching programs (ppSub and ppStack). The moments-based + // value is probably more representative of the value needed by those tools (it is also + // more stable?) + bool status = false; + float fwhmMajorMoments = psMetadataLookupF32(&status, readout->analysis, "IQ_FW1"); + float fwhmMinorMoments = psMetadataLookupF32(&status, readout->analysis, "IQ_FW1"); + + pmChip *chip = readout->parent->parent; // Parent chip + psAssert(chip, "Cell should be attached to a chip."); + psMetadataItem *item = psMetadataLookup(chip->concepts, "CHIP.SEEING"); // Item with chip + item->data.F32 = 0.5 * (fwhmMajorMoments + fwhmMinorMoments); + + psLogMsg ("psphot", PS_LOG_DETAIL, "fwhm (psf): %f,%f (moments): %f,%f", fwhmMaj, fwhmMin, fwhmMajorMoments, fwhmMinorMoments); + } + + if (psfExtra1->n) { + if (!psVectorStats (stats, psfExtra1, NULL, NULL, 0)) { + psError(PS_ERR_UNKNOWN, false, "failure to measure stats for PSF EXTRA 1"); + goto escape; + } + psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "PSF_EXT1", PS_META_REPLACE, "PSF extra param 1", stats->sampleMean); + psLogMsg ("psphot", PS_LOG_DETAIL, "PSF extra parameter 1: %f +/- %f", stats->sampleMean, stats->sampleStdev); } if (psfExtra2->n) { - if (!psVectorStats (stats, psfExtra2, NULL, NULL, 0)) { - psError(PS_ERR_UNKNOWN, false, "failure to measure stats for PSF EXTRA 2"); - goto escape; - } - psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "PSF_EXT2", PS_META_REPLACE, "PSF extra param 2", stats->sampleMean); - psLogMsg ("psphot", PS_LOG_DETAIL, "PSF extra parameter 2: %f +/- %f", stats->sampleMean, stats->sampleStdev); + if (!psVectorStats (stats, psfExtra2, NULL, NULL, 0)) { + psError(PS_ERR_UNKNOWN, false, "failure to measure stats for PSF EXTRA 2"); + goto escape; + } + psMetadataAddF32 (readout->analysis, PS_LIST_TAIL, "PSF_EXT2", PS_META_REPLACE, "PSF extra param 2", stats->sampleMean); + psLogMsg ("psphot", PS_LOG_DETAIL, "PSF extra parameter 2: %f +/- %f", stats->sampleMean, stats->sampleStdev); } Index: branches/eam_branches/ipp-20140904/psphot/src/psphotCleanup.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotCleanup.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotCleanup.c (revision 37804) @@ -25,5 +25,5 @@ // don't display memory leaks unless trace is on. // fprintf (stderr, "found %d leaks at %s\n", psMemCheckLeaks (0, NULL, NULL, false), "psphot"); - fprintf (stderr, "found %d leaks at %s\n", psMemCheckLeaks (0, NULL, stdout, false), "psphot"); + fprintf (stderr, "found %d leaks at %s\n", psMemCheckLeaks2 (0, NULL, stdout, false, 500), "psphot"); #endif return; Index: branches/eam_branches/ipp-20140904/psphot/src/psphotEfficiency.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotEfficiency.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotEfficiency.c (revision 37804) @@ -502,5 +502,5 @@ // measure the magnitudes and fluxes for the sources - if (!psphotMagnitudesReadout(config, recipe, view, readout, fakeSourcesAll, psf)) { + if (!psphotMagnitudesReadout(config, recipe, view, readout, fakeSourcesAll, psf, index)) { psError(PS_ERR_UNKNOWN, false, "Unable to measure magnitudes of fake sources."); psFree(fakeSources); Index: branches/eam_branches/ipp-20140904/psphot/src/psphotExtendedSourceFits.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotExtendedSourceFits.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotExtendedSourceFits.c (revision 37804) @@ -629,5 +629,5 @@ pmSourceCacheModel (source, maskVal); -# if (PS_TRACE_ON) +# if (0 && PS_TRACE_ON) pmModel *model = source->modelFits->data[0]; int flags = 0xffffffff; Index: branches/eam_branches/ipp-20140904/psphot/src/psphotFitSourcesLinear.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotFitSourcesLinear.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotFitSourcesLinear.c (revision 37804) @@ -263,4 +263,10 @@ } } + if (modelSum == 0 || maskedSum == 0) { + // no pixels at all drop the modelPSF (probably only an issue for external sources + // which are created from a fake peak (full force)) + psFree(source->modelPSF); + continue; + } if (modelSum < 0.5) continue; // skip sources with no model constraint (somewhat arbitrary limit) if (modelSum < 0.8) { Index: branches/eam_branches/ipp-20140904/psphot/src/psphotFullForceReadout.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotFullForceReadout.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotFullForceReadout.c (revision 37804) @@ -117,4 +117,10 @@ psphotGalaxyShape (config, view, filerule); + // measure aperture photometry corrections + if (!psphotApResid (config, view, filerule)) { + psLogMsg ("psphot", 3, "failed on psphotApResid"); + return psphotReadoutCleanup (config, view, filerule); + } + // calculate source magnitudes (psf mag and ap mag) psphotMagnitudes(config, view, filerule); Index: branches/eam_branches/ipp-20140904/psphot/src/psphotGalaxyParams.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotGalaxyParams.c (revision 37804) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotGalaxyParams.c (revision 37804) @@ -0,0 +1,637 @@ +# include "psphotInternal.h" + +#define USE_FOOTPRINTS 1 + +// measure some properties on the exended objects + +bool psphotGalaxyParams (pmConfig *config, const pmFPAview *view, const char *filerule) +{ + bool status = true; + + fprintf (stdout, "\n"); + psLogMsg ("psphot", PS_LOG_INFO, "--- psphot Galaxy Parameters ---"); + + // select the appropriate recipe information + psMetadata *recipe = psMetadataLookupPtr (&status, config->recipes, PSPHOT_RECIPE); + psAssert (recipe, "missing recipe?"); + + // add this + bool doGalaxyParams = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_GALAXY_PARAMS"); + bool doPetrosian = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_PETROSIAN"); + bool doAnnuli = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_ANNULI"); + bool doExtModels = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_FITS"); + + // All of these must be true to measure galaxy parameters + if (!doGalaxyParams || !doPetrosian || !doAnnuli || !doExtModels) { + psLogMsg ("psphot", PS_LOG_INFO, "skipping Galaxy parameter source measurements\n"); + return true; + } + + int num = psphotFileruleCount(config, filerule); + + // loop over the available readouts + for (int i = 0; i < num; i++) { + if (!psphotGalaxyParamsReadout (config, view, filerule, i, recipe)) { + psError (PSPHOT_ERR_CONFIG, false, "failed on galaxy parameter measurements for %s entry %d", filerule, i); + return false; + } + } + return true; +} + +#ifdef notdef +/*** for the moment, this test code : it is not thread safe ***/ +int Nall = 0; +int Nskip1 = 0; +int Nskip2 = 0; +int Nskip3 = 0; +int Nskip4 = 0; +int Nskip5 = 0; +int Nskip6 = 0; + +# define SKIP(VALUE) { VALUE++; continue; } + +#endif //notdef + +#define SKIP(VALUE) {continue;} + +// aperture-like measurements for extended sources +bool psphotGalaxyParamsReadout (pmConfig *config, const pmFPAview *view, const char *filerule, int index, psMetadata *recipe) { + + bool status; + int Next = 0; + int Npetro = 0; + int Nannuli = 0; + + psTimerStart ("psphot.galaxy"); + + // find the currently selected readout + pmFPAfile *file = pmFPAfileSelectSingle(config->files, filerule, index); // File of interest + psAssert (file, "missing file?"); + + pmReadout *readout = pmFPAviewThisReadout(view, file->fpa); + psAssert (readout, "missing readout?"); + + pmDetections *detections = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.DETECTIONS"); + psAssert (detections, "missing detections?"); + + psArray *sources = detections->allSources; + psAssert (sources, "missing sources?"); + + if (!sources->n) { + psLogMsg ("psphot", PS_LOG_INFO, "no sources, skipping source size"); + return true; + } + + // determine the number of allowed threads + int nThreads = psMetadataLookupS32(&status, config->arguments, "NTHREADS"); // Number of threads + if (!status) { + nThreads = 0; + } + + // get the sky noise from the background analysis; if this is missing, get the user-supplied value + float skynoise = psMetadataLookupF32 (&status, readout->analysis, "SKY_STDEV"); + if (!status) { + skynoise = psMetadataLookupF32 (&status, recipe, "SKY.NOISE"); + psWarning ("failed to get sky noise level from background analysis; defaulting to user supplied value of %f\n", skynoise); + } + + float fitNsigmaConv = psMetadataLookupF32 (&status, recipe, "EXT_FIT_NSIGMA_CONV"); // number of sigma for the convolutio + if (!status || !isfinite(fitNsigmaConv) || fitNsigmaConv <= 0) { + fitNsigmaConv = 5.0; + } + int psfSize = psMetadataLookupS32 (&status, recipe, "PCM_BOX_SIZE"); + + psImage *footprintImage = NULL; +#ifdef USE_FOOTPRINTS + footprintImage = psImageAlloc(readout->image->numCols, readout->image->numRows, PS_TYPE_S32); + psImageInit(footprintImage, 0); + pmSetFootprintArrayIDsForImage(footprintImage, detections->footprints, true); +#endif + + + // source analysis is done in S/N order (brightest first) + sources = psArraySort (sources, pmSourceSortByFlux); + + // option to limit analysis to a specific region + char *region = psMetadataLookupStr (&status, recipe, "ANALYSIS_REGION"); + psRegion *AnalysisRegion = psRegionAlloc(0,0,0,0); + *AnalysisRegion = psRegionForImage (readout->image, psRegionFromString (region)); + if (psRegionIsNaN (*AnalysisRegion)) psAbort("analysis region mis-defined"); + + // choose Cx, Cy (see psphotThreadTools.c for overview of the concepts) + int Cx = 1, Cy = 1; + psphotChooseCellSizes (&Cx, &Cy, readout, nThreads); + + psArray *cellGroups = psphotAssignSources (Cx, Cy, sources); + + for (int i = 0; i < cellGroups->n; i++) { + + psArray *cells = cellGroups->data[i]; + + for (int j = 0; j < cells->n; j++) { + + // allocate a job -- if threads are not defined, this just runs the job + psThreadJob *job = psThreadJobAlloc ("PSPHOT_GALAXY_PARAMS"); + + psArrayAdd(job->args, 1, readout); + psArrayAdd(job->args, 1, cells->data[j]); // sources + psArrayAdd(job->args, 1, AnalysisRegion); + psArrayAdd(job->args, 1, recipe); + + PS_ARRAY_ADD_SCALAR(job->args, skynoise, PS_TYPE_F32); + PS_ARRAY_ADD_SCALAR(job->args, fitNsigmaConv, PS_TYPE_F32); + PS_ARRAY_ADD_SCALAR(job->args, psfSize, PS_TYPE_F32); + psArrayAdd(job->args, 1, footprintImage); + + PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Next + PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Npetro + PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nannuli + + +// set this to 0 to run without threading +# if (0) + if (!psThreadJobAddPending(job)) { + psError(PS_ERR_UNKNOWN, false, "Unable to guess model."); + psFree(AnalysisRegion); + return false; + } +# else + if (!psphotGalaxyParams_Threaded(job)) { + psError(PS_ERR_UNKNOWN, false, "Unable to guess model."); + psFree(AnalysisRegion); + return false; + } + psScalar *scalar = NULL; + scalar = job->args->data[8]; + Next += scalar->data.S32; + scalar = job->args->data[9]; + Npetro += scalar->data.S32; + scalar = job->args->data[10]; + Nannuli += scalar->data.S32; + psFree(job); +# endif + } + + // wait for the threads to finish and manage results + if (!psThreadPoolWait (false, true)) { + psError(PS_ERR_UNKNOWN, false, "Unable to guess model."); + psFree(AnalysisRegion); + return false; + } + + // we have only supplied one type of job, so we can assume the types here + psThreadJob *job = NULL; + while ((job = psThreadJobGetDone()) != NULL) { + if (job->args->n < 1) { + fprintf (stderr, "error with job\n"); + } else { + psScalar *scalar = NULL; + scalar = job->args->data[8]; + Next += scalar->data.S32; + scalar = job->args->data[9]; + Npetro += scalar->data.S32; + scalar = job->args->data[10]; + Nannuli += scalar->data.S32; + } + psFree(job); + } + } + psFree (cellGroups); + psFree(AnalysisRegion); + psFree(footprintImage); + + psLogMsg ("psphot", PS_LOG_WARN, "galaxy parameter analysis: %f sec for %d objects\n", psTimerMark ("psphot.galaxy"), Next); + psLogMsg ("psphot", PS_LOG_INFO, " %d petrosian\n", Npetro); + psLogMsg ("psphot", PS_LOG_INFO, " %d annuli\n", Nannuli); + +// # if (PS_TRACE_ON) +# if (0) + fprintf (stderr, "ext analysis skipped @ 1 : %d\n", Nskip1); + fprintf (stderr, "ext analysis skipped @ 2 : %d\n", Nskip2); + fprintf (stderr, "ext analysis skipped @ 3 : %d\n", Nskip3); + fprintf (stderr, "ext analysis skipped @ 4 : %d\n", Nskip4); + fprintf (stderr, "ext analysis skipped @ 5 : %d\n", Nskip5); + fprintf (stderr, "ext analysis skipped @ 6 : %d\n", Nskip6); +#endif + + // psphotRadialProfileShowSkips (); + + // psphotVisualShowResidualImage (readout, false); + + // psphotVisualShowPetrosians (sources); + + return true; +} + +bool psphotGalaxyParams_Threaded (psThreadJob *job) { + + // bool status; + + int Next = 0; + int Npetro = 0; + int Nannuli = 0; + int Ngood = 0; + + // arguments: readout, sources, models, region, psfSize, maskVal, markVal + pmReadout *readout = job->args->data[0]; + psArray *sources = job->args->data[1]; + psRegion *region = job->args->data[2]; + psMetadata *recipe = job->args->data[3]; + +#ifndef USE_FOOTPRINTS + float skynoise = PS_SCALAR_VALUE(job->args->data[4],F32); +#endif + float fitNsigmaConv = PS_SCALAR_VALUE(job->args->data[5],F32); + float psfSize = PS_SCALAR_VALUE(job->args->data[6],F32); + psImage *footprintImage = job->args->data[7]; + + bool status; +// bool doPetrosian = psMetadataLookupBool (&status, recipe, "EXTENDED_SOURCE_PETROSIAN"); + + // user-defined masks to test for good/bad pixels (build from recipe list if not yet set) + psImageMaskType maskVal = psMetadataLookupImageMask(&status, recipe, "MASK.PSPHOT"); // Mask value for bad pixels + assert (maskVal); + + // don't try and use bad models + pmModelStatus badModel = PM_MODEL_STATUS_NONE; + badModel |= PM_MODEL_STATUS_NONCONVERGE; + badModel |= PM_MODEL_STATUS_BADARGS; + badModel |= PM_MODEL_STATUS_OFFIMAGE; + badModel |= PM_MODEL_STATUS_NAN_CHISQ; + badModel |= PM_MODEL_SERSIC_PCM_FAIL_GUESS; + badModel |= PM_MODEL_SERSIC_PCM_FAIL_GRID; + badModel |= PM_MODEL_PCM_FAIL_GUESS; + + pmModelType sersicModelType = pmModelClassGetType("PS_MODEL_SERSIC"); + + psVector *bgPixels = NULL; + psVector *bgPixelsFlip = NULL; +#define dontMAKE_MARK_IMAGE 1 +#ifdef MAKE_MARK_IMAGE + psImage *markImage = NULL; +#endif + + #define MARK_SOURCE 0 + #define MARK_MASK 1 + #define MARK_BG 2 + psRandom *rng = psRandomAlloc(PS_RANDOM_TAUS); + +#define dontDUMP_IMAGES 1 +#ifdef DUMP_IMAGES + char *filebase = "/data/ippc18.0/home/bills/debug/psphot/images"; + char filename[80]; + psMetadata *dummyHeader = psMetadataAlloc(); +#endif + + // choose the sources of interest + for (int i = 0; i < sources->n; i++) { + + pmSource *source = sources->data[i]; + + // Nall ++; + + // rules for measuring petrosian parameters for specific objects are set in + // psphotChooseAnalysisOptions.c + if (!(source->tmpFlags & PM_SOURCE_TMPF_PETRO)) SKIP (Nskip1); + + // limit selection by analysis region (XXX move this into psphotChooseAnalysisOption?) + if (source->peak->x < region->x0) SKIP (Nskip2); + if (source->peak->y < region->y0) SKIP (Nskip3); + if (source->peak->x > region->x1) SKIP (Nskip4); + if (source->peak->y > region->y1) SKIP (Nskip5); + + // Should we skip objects that do not have an extended model? + if (!source->extpars) SKIP (Nskip6); + if (!source->modelFits) SKIP (Nskip7); + + // Flux cut to reduce sources to those that match SDSS + // if (source->moments->KronFlux < 2e5) continue; + + pmModel *sersicModel = NULL; + for (int j = 0; j < source->modelFits->n; j++) { + pmModel *model = source->modelFits->data[j]; + if (model->type == sersicModelType) { + sersicModel = model; + break; + } + } + if (!sersicModel) SKIP (Nskip8); + if (sersicModel->flags & badModel) SKIP (Nskip8); + Next ++; + + // if sersicModel is not modelEXT make it so + pmSourceMode saveMode = source->mode; + pmSourceType saveType = source->type; + bool isExtended = saveType == PM_SOURCE_TYPE_EXTENDED; + pmModel *extModel = source->modelEXT; + + if (source->tmpFlags & PM_SOURCE_TMPF_SUBTRACTED) { + pmSourceAdd (source, PM_MODEL_OP_FULL, maskVal); + } + + if (!isExtended || extModel != sersicModel) { + // XXX: I'm not doing this right. For now just look at objects for which + source->modelEXT = psMemIncrRefCounter (sersicModel); + source->type = PM_SOURCE_TYPE_EXTENDED; + source->mode |= PM_SOURCE_MODE_EXTMODEL; + // do we need this? + source->mode |= PM_SOURCE_MODE_NONLINEAR_FIT; +#ifdef notdef +#endif + } + + + psF32 *PAR = sersicModel->params->data.F32; + psF32 Xo = PAR[PM_PAR_XPOS]; + psF32 Yo = PAR[PM_PAR_YPOS]; + psF32 Io = PAR[PM_PAR_I0]; + // For now: set the half light radius to the major axis of the model + psEllipseAxes axes = pmPSF_ModelToAxes (PAR, sersicModel->class->useReff); + source->extpars->ghalfLightRadius = axes.major; + // and we will integrate out to 2 * half light radius + psF32 R = 2 * axes.major; + psF32 R2 = PS_SQR(R); + + float radius = source->peak->xf - source->pixels->col0; + radius = PS_MAX (radius, source->peak->yf - source->pixels->row0); + radius = PS_MAX (radius, source->pixels->numRows - source->peak->yf + source->pixels->row0); + radius = PS_MAX (radius, source->pixels->numCols - source->peak->xf + source->pixels->col0); + // Make image at least 10 * half light radius + radius = PS_MAX (radius, 5*axes.major); + pmSourceRedefinePixels (source, readout, source->peak->xf, source->peak->yf, 2*radius); + + // RedefinePixels frees modelFlux + if (source->modelEXT->isPCM) { + pmPCMCacheModel (source, maskVal, psfSize, fitNsigmaConv); + } else { + pmSourceCacheModel (source, maskVal); + } + +#ifdef DUMP_IMAGES + sprintf(filename, "%s/%04d-pix.fits", filebase, source->id); + + psphotSaveImage(dummyHeader, source->pixels, filename); +#endif + + + // Create the markImage and accumlate the Flux + + // center of mass in subimage. Note: the calculation below uses pixel index, so we correct + // xCM, yCM from pixel coords to pixel index here. + psF32 xCM = Xo - 0.5 - source->pixels->col0; // coord of peak in subimage + psF32 yCM = Yo - 0.5 - source->pixels->row0; // coord of peak in subimage + + psF32 **vPix = source->pixels->data.F32; + psImageMaskType **vMsk = (source->maskObj == NULL) ? NULL : source->maskObj->data.PS_TYPE_IMAGE_MASK_DATA; + + bgPixels = psVectorRecycle(bgPixels, source->pixels->numCols*source->pixels->numRows, PS_TYPE_F32); + bgPixels->n = 0; + bgPixelsFlip = psVectorRecycle(bgPixelsFlip, source->pixels->numCols*source->pixels->numRows, PS_TYPE_F32); + bgPixelsFlip->n = 0; + +#ifdef MAKE_MARK_IMAGE + // The markImage may no longer be needed + markImage = psImageRecycle(markImage, source->pixels->numCols, source->pixels->numRows, PS_TYPE_U8); + + // assume all pixels are masked to start + psImageInit(markImage, MARK_MASK); + psU8 **vMark = markImage->data.U8; + #define SET_MARK(_val) vMark[row][col] = _val +#else + #define SET_MARK(_val) +#endif + +#ifdef USE_FOOTPRINTS + psS32 **vFootprint = footprintImage->data.S32; +#else + // XXX use recipe value instead of hard coded 1.5 + psF32 threshold = skynoise * 3; +#endif + + psF32 sumI = 0; + psS32 row0 = source->pixels->row0; + psS32 col0 = source->pixels->col0; + for (psS32 row = 0; row < source->pixels->numRows ; row++) { + psF32 yDiff = row - yCM; + // y coordinate of mirror pixel + int yFlip = yCM - yDiff; + if (yFlip < 0) continue; + if (yFlip >= source->pixels->numRows) continue; + psS32 yImage = row + row0; + for (psS32 col = 0; col < source->pixels->numCols ; col++) { + // check mask and value for this pixel + if (vMsk && (vMsk[row][col] & maskVal)) continue; + psF32 pixel = vPix[row][col]; + if (isnan(pixel)) continue; + // valid pixel + + psS32 xImage = col + col0; + + psF32 xDiff = col - xCM; + + psF32 r2 = PS_SQR(xDiff) + PS_SQR(yDiff); + if (r2 <= R2) { + sumI += pixel; + } + +#ifdef USE_FOOTPRINTS + if (vFootprint[yImage][xImage]) { +#else + if (fabs(pixel) >= threshold) { +#endif + SET_MARK(MARK_SOURCE); + continue; + } + // we have a background pixel. + // Check whether it's mirror is a background pixel as well. + + // x coordinate of mirror pixel + int xFlip = xCM - xDiff; + if (xFlip < 0) continue; + if (xFlip >= source->pixels->numCols) continue; + // check mask and value for mirror pixel + if (vMsk && (vMsk[yFlip][xFlip] & maskVal)) continue; + psF32 mirrorPixel = vPix[yFlip][xFlip]; + if (isnan(mirrorPixel)) continue; + +#ifdef USE_FOOTPRINTS + if (vFootprint[yFlip + row0][xFlip + col0] == 0) { +#else + if (fabs(mirrorPixel) < threshold) { +#endif + // This is a good background pixel with a good mirror Save the values + SET_MARK(MARK_BG); + psVectorAppend(bgPixels, pixel); + psVectorAppend(bgPixelsFlip, mirrorPixel); + } + } + } +#ifdef DUMP_IMAGES + sprintf(filename, "%s/%04d-mrk.fits", filebase, source->id); + + psphotSaveImage(dummyHeader, markImage, filename); +#endif + + // subtract model to create the residual image + pmSourceSub (source, PM_MODEL_OP_FULL, maskVal); + if (!(source->tmpFlags & PM_SOURCE_TMPF_SUBTRACTED)) { + // Hmm something has gone wrong + fprintf(stderr, "failed to subtract sersic model for source %d\n", source->id); + goto restoreModel; + } +#ifdef DUMP_IMAGES + sprintf(filename, "%s/%04d-res.fits", filebase, source->id); + + psphotSaveImage(dummyHeader, source->pixels, filename); +#endif + + // don't know if these pointers could change in the background subtraction so grab them again + vPix = source->pixels->data.F32; + vMsk = (source->maskObj == NULL) ? NULL : source->maskObj->data.PS_TYPE_IMAGE_MASK_DATA; + psF32 **vVar = source->variance->data.F32; + psF32 **vModel = source->modelFlux->data.F32; + // Accumulate Rt and Ra + psF32 sumRt = 0; + psF32 sumRa = 0; + int numPixels = 0; + psF32 sumB = 0; + psF32 sumModel = 0; + int numB = 0; + for (psS32 row = 0; row < source->pixels->numRows ; row++) { + psF32 yDiff = row - yCM; + if (fabs(yDiff) > R) continue; + // y coordinate of mirror pixel + int yFlip = yCM - yDiff; + if (yFlip < 0) continue; + if (yFlip >= source->pixels->numRows) continue; + for (psS32 col = 0; col < source->pixels->numCols ; col++) { + // check mask and value for this pixel + if (vMsk && (vMsk[row][col] & maskVal)) continue; + psF32 pixel = vPix[row][col]; + if (isnan(pixel)) continue; + + psF32 xDiff = col - xCM; + psF32 r2 = PS_SQR(xDiff) + PS_SQR(yDiff); + if (r2 > R2) continue; + + // Measurement of bumpiness parameter excludes the pixels + // within 2 pixels of the center + if (r2 > 4) { + psF32 variance = vVar[row][col]; + if (isfinite(variance)) { + psF32 modelPixel = vModel[row][col]; + if (!isfinite(modelPixel)) { + // XXX: Can this happen in a good model? + fprintf(stderr, "non-finite modelPixel %4d (%4d, %4d)\n", source->id, col, row); + } + numB++; + sumB += PS_SQR(pixel) - variance; + sumModel += Io * modelPixel; + } + } + + // x coordinate of mirror pixel + int xFlip = xCM - xDiff; + if (xFlip < 0) continue; + if (xFlip >= source->pixels->numCols) continue; + // check mask and value for mirror pixel + if (vMsk && (vMsk[yFlip][xFlip] & maskVal)) continue; + psF32 mirrorPixel = vPix[yFlip][xFlip]; + if (isnan(mirrorPixel)) continue; + + sumRt += fabs(pixel + mirrorPixel); + sumRa += fabs(pixel - mirrorPixel); + ++numPixels; + } + } + + if (numB > 0) { + source->extpars->gbumpy = sqrt(sumB/numB) / (sumModel / numB); + } + + if (bgPixels->n < numPixels) { + fprintf(stderr, "could not find enough background pixels for source %d found %ld need %d\n", source->id, bgPixels->n, numPixels); + goto restoreModel; + } + Ngood++; + + // and calculate Bt and Ba + psF32 sumBt = 0; + psF32 sumBa = 0; + psF32 *vB = bgPixels->data.F32; + psF32 *vBFlip = bgPixelsFlip->data.F32; + long length = bgPixels->n; + // Randomly select numPixels background pixels + for (int n = 0; n < numPixels; n++) { + double rand = psRandomUniform(rng); + long index = floor(rand * length); + if (index >= length) { + fprintf(stderr, "random index %ld is larger than vector length %ld\n", + index, length); + continue; + } + sumBt += fabs(vB[index] + vBFlip[index]); + sumBa += fabs(vB[index] - vBFlip[index]); + } + + // finally save the results + source->extpars->gRT = 0.5 * (sumRt - sumBt)/sumI; + source->extpars->gRA = 0.5 * (sumRa - sumBa)/sumI; + source->extpars->gS2 = source->extpars->gRT + source->extpars->gRA; + +#ifdef PRINT_STUFF + if (Ngood % 40 == 1) { + printf(" id sumI rHL rad sumI gRT gRA sumRt sumBt numPixels numBGPixels\n"); + } + printf("%5d %5.2e %3.1f %10.1f %6.3f %6.3f %5.2e %5.2e %10d %6ld\n", + source->id, sumI, source->extpars->ghalfLightRadius, radius, source->extpars->gRT, source->extpars->gRA, sumRt, sumBt, numPixels, bgPixels->n); +#endif + + // Now put things back the way that they were. +restoreModel: + + if ((extModel != sersicModel) || !isExtended) { + // add the sersic model back in + if (source->tmpFlags & PM_SOURCE_TMPF_SUBTRACTED) { + pmSourceAdd (source, PM_MODEL_OP_FULL, maskVal); + } + + // replace the extended model + if (extModel != sersicModel) { + psMemDecrRefCounter (sersicModel); + source->modelEXT = extModel; + } + source->mode = saveMode; + source->type = saveType; + pmSourceSub (source, PM_MODEL_OP_FULL, maskVal); + } + } + +#ifdef MAKE_MARK_IMAGE + psFree(markImage); +#endif +#ifdef DUMP_IMAGES + psFree(dummyHeader); + // psphotSaveImage raises an error + psErrorClear(); +#endif + psFree(bgPixels) + psFree(bgPixelsFlip) + psFree(rng); + + psScalar *scalar = NULL; + + // change the value of a scalar on the array (wrap this and put it in psArray.h) + scalar = job->args->data[8]; + scalar->data.S32 = Next; + + scalar = job->args->data[9]; + scalar->data.S32 = Npetro; + + scalar = job->args->data[10]; + scalar->data.S32 = Nannuli; + + return true; +} + Index: branches/eam_branches/ipp-20140904/psphot/src/psphotMagnitudes.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotMagnitudes.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotMagnitudes.c (revision 37804) @@ -38,5 +38,5 @@ psAssert (psf, "missing psf?"); - if (!psphotMagnitudesReadout (config, recipe, view, readout, sources, psf)) { + if (!psphotMagnitudesReadout (config, recipe, view, readout, sources, psf, i)) { psError (PSPHOT_ERR_CONFIG, false, "failed to measure magnitudes for %s entry %d", filerule, i); return false; @@ -46,5 +46,5 @@ } -bool psphotMagnitudesReadout(pmConfig *config, psMetadata *recipe, const pmFPAview *view, pmReadout *readout, psArray *sources, pmPSF *psf) { +bool psphotMagnitudesReadout(pmConfig *config, psMetadata *recipe, const pmFPAview *view, pmReadout *readout, psArray *sources, pmPSF *psf, int index) { bool status = false; @@ -65,6 +65,6 @@ // if backModel or backStdev are missing, the values of sky and/or skyErr will be set to NAN - pmReadout *backModel = psphotSelectBackground (config, view); - pmReadout *backStdev = psphotSelectBackgroundStdev (config, view); + pmReadout *backModel = psphotSelectBackground (config, view, index); + pmReadout *backStdev = psphotSelectBackgroundStdev (config, view, index); // bit-masks to test for good/bad pixels Index: branches/eam_branches/ipp-20140904/psphot/src/psphotOutput.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotOutput.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotOutput.c (revision 37804) @@ -29,10 +29,7 @@ } -pmReadout *psphotSelectBackground (pmConfig *config, const pmFPAview *view) { - - bool status; - - - pmFPAfile *file = psMetadataLookupPtr (&status, config->files, psphotGetFilerule("PSPHOT.BACKMDL")); +pmReadout *psphotSelectBackground (pmConfig *config, const pmFPAview *view, int index) { + + pmFPAfile *file = pmFPAfileSelectSingle (config->files, psphotGetFilerule("PSPHOT.BACKMDL"), index); if (!file) return NULL; @@ -41,9 +38,7 @@ } -pmReadout *psphotSelectBackgroundStdev (pmConfig *config, const pmFPAview *view) { - - bool status; - - pmFPAfile *file = psMetadataLookupPtr (&status, config->files, psphotGetFilerule("PSPHOT.BACKMDL.STDEV")); +pmReadout *psphotSelectBackgroundStdev (pmConfig *config, const pmFPAview *view, int index) { + + pmFPAfile *file = pmFPAfileSelectSingle (config->files, psphotGetFilerule("PSPHOT.BACKMDL.STDEV"), index); if (!file) return NULL; Index: branches/eam_branches/ipp-20140904/psphot/src/psphotPetrosianAnalysis.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotPetrosianAnalysis.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotPetrosianAnalysis.c (revision 37804) @@ -12,11 +12,4 @@ // XXX temporary user-supplied systematic sky noise measurement (derive from background model) float skynoise = psMetadataLookupF32 (&status, recipe, "SKY.NOISE"); - -# if (0) - // if backModel or backStdev are missing, the values of sky and/or skyErr will be set to NAN - // XXX use this to set skynoise - pmReadout *backModel = psphotSelectBackground (config, view); - pmReadout *backStdev = psphotSelectBackgroundStdev (config, view); -# endif // S/N limit to perform full non-linear fits Index: branches/eam_branches/ipp-20140904/psphot/src/psphotReadout.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotReadout.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotReadout.c (revision 37804) @@ -325,4 +325,6 @@ psphotExtendedSourceAnalysis (config, view, filerule); // pass 1 (detections->allSources) psphotExtendedSourceFits (config, view, filerule); // pass 1 (detections->allSources) + // measure some parameters for galaxy science + psphotGalaxyParams (config, view, filerule); psphotRadialApertures(config, view, filerule, 0); Index: branches/eam_branches/ipp-20140904/psphot/src/psphotReadoutCleanup.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotReadoutCleanup.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotReadoutCleanup.c (revision 37804) @@ -57,5 +57,5 @@ // use the psf-model to measure FWHM stats if (psf) { - if (!psphotPSFstats (readout, psf)) { + if (!psphotPSFstatsSources (readout, sources, psf)) { psError(PSPHOT_ERR_PROG, false, "Failed to measure PSF shape parameters"); return false; Index: branches/eam_branches/ipp-20140904/psphot/src/psphotSetNFrames.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotSetNFrames.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotSetNFrames.c (revision 37804) @@ -45,8 +45,14 @@ pmDetections *detections = psMetadataLookupPtr (&status, readout->analysis, "PSPHOT.DETECTIONS"); - psAssert (detections, "missing detections?"); - psAssert (detections->allSources, "missing sources?"); + if (!detections) { + psLogMsg ("psphot", PS_LOG_INFO, "No detections image input %d", index); + return true; + } psArray *sources = detections->allSources; + if (!sources) { + psLogMsg ("psphot", PS_LOG_INFO, "No sources for input %d", index); + return true; + } int col0 = expnum->col0; Index: branches/eam_branches/ipp-20140904/psphot/src/psphotSetThreads.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotSetThreads.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotSetThreads.c (revision 37804) @@ -92,4 +92,9 @@ psFree(task); + task = psThreadTaskAlloc("PSPHOT_GALAXY_PARAMS", 11); + task->function = &psphotGalaxyParams_Threaded; + psThreadTaskAdd(task); + psFree(task); + return true; } Index: branches/eam_branches/ipp-20140904/psphot/src/psphotSourceStats.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotSourceStats.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotSourceStats.c (revision 37804) @@ -185,18 +185,18 @@ psThreadJob *job = psThreadJobAlloc ("PSPHOT_SOURCE_STATS"); - psArrayAdd(job->args, 1, cells->data[j]); // sources - - PS_ARRAY_ADD_SCALAR(job->args, INNER, PS_TYPE_F32); - PS_ARRAY_ADD_SCALAR(job->args, MIN_SN, PS_TYPE_F32); - PS_ARRAY_ADD_SCALAR(job->args, RADIUS, PS_TYPE_F32); - PS_ARRAY_ADD_SCALAR(job->args, SIGMA, PS_TYPE_F32); - PS_ARRAY_ADD_SCALAR(job->args, MIN_KRON_RADIUS, PS_TYPE_F32); - - PS_ARRAY_ADD_SCALAR(job->args, maskVal, PS_TYPE_IMAGE_MASK); - PS_ARRAY_ADD_SCALAR(job->args, markVal, PS_TYPE_IMAGE_MASK); - - PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nmoments - PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nfail - PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // this is used as a return value for Nfaint + psArrayAdd(job->args, 1, cells->data[j]); // 0 : sources + + PS_ARRAY_ADD_SCALAR(job->args, INNER, PS_TYPE_F32); // 1 + PS_ARRAY_ADD_SCALAR(job->args, MIN_SN, PS_TYPE_F32); // 2 + PS_ARRAY_ADD_SCALAR(job->args, RADIUS, PS_TYPE_F32); // 3 + PS_ARRAY_ADD_SCALAR(job->args, SIGMA, PS_TYPE_F32); // 4 + PS_ARRAY_ADD_SCALAR(job->args, MIN_KRON_RADIUS, PS_TYPE_F32); // 5 + + PS_ARRAY_ADD_SCALAR(job->args, maskVal, PS_TYPE_IMAGE_MASK); // 6 + PS_ARRAY_ADD_SCALAR(job->args, markVal, PS_TYPE_IMAGE_MASK); // 7 + + PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // 8 : this is used as a return value for Nmoments + PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // 9 : this is used as a return value for Nfail + PS_ARRAY_ADD_SCALAR(job->args, 0, PS_TYPE_S32); // 10 : this is used as a return value for Nfaint if (!psThreadJobAddPending(job)) { @@ -316,4 +316,5 @@ // XXX: this must match the above + psAbort ("this code is broken: must match psphotSourceStatsReadout above"); psArrayAdd(job->args, 1, cells->data[j]); // sources psArrayAdd(job->args, 1, recipe); @@ -472,11 +473,11 @@ // change the value of a scalar on the array (wrap this and put it in psArray.h) - scalar = job->args->data[7]; + scalar = job->args->data[8]; scalar->data.S32 = Nmoments; - scalar = job->args->data[8]; + scalar = job->args->data[9]; scalar->data.S32 = Nfail; - scalar = job->args->data[9]; + scalar = job->args->data[10]; scalar->data.S32 = Nfaint; Index: branches/eam_branches/ipp-20140904/psphot/src/psphotStackParseCamera.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotStackParseCamera.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotStackParseCamera.c (revision 37804) @@ -262,4 +262,5 @@ } +#ifdef MAKE_CHISQ_IMAGE // generate an pmFPAimage for the chisqImage { @@ -293,4 +294,7 @@ chisqVariance->save = saveChisq; } +#else + (void) saveChisq; +#endif // MAKE_CHISQ_IMAGE psTrace("psphot", 1, "Done with psphotStackParseCamera...\n"); Index: branches/eam_branches/ipp-20140904/psphot/src/psphotStackReadout.c =================================================================== --- branches/eam_branches/ipp-20140904/psphot/src/psphotStackReadout.c (revision 37803) +++ branches/eam_branches/ipp-20140904/psphot/src/psphotStackReadout.c (revision 37804) @@ -341,4 +341,7 @@ // measure non-linear extended source models (exponential, deVaucouleur, Sersic) (sources sorted by S/N) psphotExtendedSourceFits (config, view, STACK_RAW); // pass 1 (detections->allSources) + + // measure some parameters for galxy science + psphotGalaxyParams (config, view, STACK_RAW); // pass 1 (detections->allSources) // create source children for the OUT filerule (for radial aperture photometry and output) @@ -354,4 +357,7 @@ return psphotReadoutCleanup (config, view, STACK_RAW); } + // These arrays are no longer needed. The inputs' sources arrays have references to the data + psFree (objects); + psFree (objectsOut); bool radial_apertures = psMetadataLookupBool(NULL, recipe, "RADIAL_APERTURES"); @@ -410,6 +416,4 @@ // measure aperture photometry corrections if (!psphotApResid (config, view, STACK_RAW)) { - psFree (objects); - psFree (objectsOut); psLogMsg ("psphot", 3, "failed on psphotApResid"); return psphotReadoutCleanup (config, view, STACK_RAW); @@ -432,5 +436,5 @@ logMemStats("final"); -#if (1) +#if (0) psphotSourceMemory(config, view, STACK_RAW); psphotSourceMemory(config, view, STACK_OUT); @@ -448,7 +452,4 @@ psphotStackRemoveChisqFromInputs(config, STACK_RAW); #endif - - psFree (objects); - psFree (objectsOut); // create the exported-metadata and free local data Index: branches/eam_branches/ipp-20140904/pstamp/scripts/pstamp_checkdependent.pl =================================================================== --- branches/eam_branches/ipp-20140904/pstamp/scripts/pstamp_checkdependent.pl (revision 37803) +++ branches/eam_branches/ipp-20140904/pstamp/scripts/pstamp_checkdependent.pl (revision 37804) @@ -24,4 +24,7 @@ my $ra_max_3PI_hours = 22.5; # disable updates for skycells with ra greater than this (hours) my $ra_max_3PI = $ra_max_3PI_hours * 15 * 3.14149 / 180.; # rawExp.ra is in radians + +# verbose flag for parse commands +my $parse_verbose = 0; # XXX: put this in a module somewhere @@ -49,4 +52,5 @@ 'dbserver=s' => \$ps_dbserver, # postage stamp server dbserver 'verbose' => \$verbose, + 'parse-verbose' => \$parse_verbose, 'save-temps' => \$save_temps, 'no-update' => \$no_update, @@ -78,4 +82,6 @@ } +print "\n==== Starting dependency checking for $dep_id $stage $stage_id $component at " . (scalar localtime) . ". ====\n"; + if ($label) { # rlabel is deprecated. Use one based on the supplied label parameter which is the current label @@ -84,4 +90,5 @@ # (and formerly in pstampparse.pl) is not particularly clean but it's simple. my $new_rlabel = 'ps_ud_' . $label; + $rlabel = "undefined" if !defined $rlabel; if ($new_rlabel ne $rlabel) { print "Notice: using $new_rlabel instead of $rlabel for update label.\n"; @@ -120,5 +127,5 @@ } -my $it = runToolAndParseExpectOne($cmd, $verbose); +my $it = runToolAndParseExpectOne($cmd, $parse_verbose); my_die("no components found", $PS_EXIT_PROG_ERROR) if ( !$it); @@ -136,4 +143,5 @@ } + # magic is no longer rquired $need_magic = 0; @@ -142,4 +150,5 @@ if ((($it->{state} eq 'full') or ($it->{state} eq 'update')) and ($it->{data_state} eq 'full')) { + print "\nDependency Satisfied for $stage $stage_id $component\n"; # This Dependency is satisfied. All done. Release the pstampJobs # @@ -178,5 +187,5 @@ } - # assume the fault is transient. + # hope the fault is transient. # fault the dependent to give the fault a chance to correct itself my_die("Component faulted on update dep_id: $dep_id", $PS_EXIT_SYS_ERROR); @@ -279,5 +288,5 @@ ($fault eq $PSTAMP_GONE)) { - print "chipRun state is $chip->{chip_id} has state: $state data_state: $data_state fault: $fault cannot update\n"; + print "chipRun state is $chip->{chip_id} has state: $state class_id $chip->{class_id} data_state: $data_state fault: $fault cannot update.\n"; my $error_code; if (($state eq 'error_cleaned') or ($data_state ='error_cleaned')) { @@ -291,10 +300,10 @@ } elsif ($chip->{state} eq 'goto_cleaned') { - # we don't need to check for full data_state - # if ($chip->{data_state} eq 'full') { { + print "dependent chip is in state goto_cleaned, changing state to update.\n"; + # cleanup must not be running. Set state to update. If this chip is not 'full' it will be # set to be updated the next time this script is invoked - # XXX: In the very rare case that this chip run is actually in the running cleanup pantasks + # XXX: In the very unlikely case that this chip run is actually in the running cleanup pantasks # queue things may get confused but we can live with that my $command = "$chiptool -updaterun -set_state update -chip_id $chip_id"; @@ -313,4 +322,6 @@ } } elsif (($chip->{data_state} ne 'update') and ($chip->{data_state} ne 'full')) { + + print "Setting chip imfile $chip_id $chip->{class_id} to be updated.\n"; # chiptool does more state checking to insure this isn't done prematurely. @@ -343,4 +354,5 @@ # we've had a number of runs in this limbo state + print "Changing chip run $chip_id to state update for $chip->{class_id} which is in data_state update.\n"; my $command = "$chiptool -updaterun -set_state update -chip_id $chip_id"; $command .= " -set_label $rlabel" if $rlabel; @@ -364,5 +376,5 @@ return $PSTAMP_GONE; } elsif ($fault_count > $max_fault_count) { - print "$stage $stage_id has faulted $fault_count times. Giving up\n"; + print "Dependency for $stage $stage_id has been found in fault state $fault_count times. Giving up. Fault: $fault\n"; # XXX: stop faulting chips for now # $fault = $PSTAMP_GONE; @@ -421,5 +433,5 @@ if (($state =~ /error/) or ($state =~ /purged/) or ($state =~ /scrubbed/) or ($state eq 'drop') or ($data_state =~ /error/) or ($data_state =~ /purged/) or ($data_state =~ /scrubbed/) or ($data_state eq 'drop')) { - print STDERR "warpRun $warp_id $skycell_id has state $state $data_state faulting jobs\n"; + print STDERR "warpRun $warp_id $skycell_id has state $state $data_state faulting jobs.\n"; my $error_code; if (($state eq 'error_cleaned') or ($data_state eq 'error_cleaned')) { @@ -430,12 +442,14 @@ return $error_code } - if (($data_state eq 'full') and ($metadata->{fault})) { + if ((($data_state eq 'full') or ($data_state eq 'update')) and ($metadata->{fault})) { # fault dependent. my $fault = $metadata->{fault}; - print STDERR "warp $warp_id $skycell_id faulted: $fault"; + print STDERR "warp $warp_id $skycell_id faulted with code: $fault."; + # XXX: TODO maybe: If ($fault != 2 and $fault != $PSTAMP_GONE) revert the skycell. The warp task only reverts fault 2 + # and sometimes warps will fault due to nfs problems that end up with fault = 4;. return $fault; } if ($metadata->{quality} ne 0) { - print STDERR "warp $warp_id $skycell_id bad quality on update: $metadata->{quality}\n"; + print STDERR "warp $warp_id $skycell_id bad quality on update: $metadata->{quality}.\n"; return $PSTAMP_GONE; } @@ -448,4 +462,5 @@ # It's "safe" to set the state to update. If the skycell is not full it will be set to update # the next time this task runs + print "Change state of warp $warp_id to update.\n"; my $command = "$warptool -updaterun -set_state update -warp_id $warp_id"; $command .= " -set_label $rlabel" if $rlabel; @@ -465,10 +480,10 @@ # get the list of input chips for this skycell my $command = "$warptool -scmap -warp_id $warp_id -skycell_id $skycell_id"; - my $data = runToolAndParse($command, $verbose); + my $data = runToolAndParse($command, $parse_verbose); if (!$data or scalar @$data == 0) { # This happens if the chipProcessedImfile disappears which happened when earlier # versions of chiptool -revertprocessedimfile didn't check the chipRun.state before # deleting the row. - print STDERR "failed to find warpSkyCelllMap for warpRun $warp_id skycell_id $skycell_id"; + print STDERR "failed to find warpSkyCellMap for warpRun $warp_id skycell_id $skycell_id"; return $PSTAMP_GONE; } @@ -478,4 +493,5 @@ my $chip_id; foreach my $chip (@$data) { + # XXX: change tools to include cam_state (camRun.state). If it is not full the warp updates will never run. my $cam_fault = $chip->{cam_fault}; if (defined $cam_fault and $cam_fault > 0) { @@ -488,7 +504,7 @@ $chip->{decl} = $metadata->{decl}; +if (0) { # if chip has been magicked before require it to be magicked again # because the warp pending query requires it. -if (0) { if ($chip->{magicked} < 0) { print "Input has been destreaked so we must destreak before warping\n"; @@ -514,4 +530,5 @@ # the reason we defer setting the warp to update is so that we can handle error conditions at previous # stages more easily. + print "Chips are ready for warp $warp_id $skycell->{skycell_id}. Setting skycell to be updated.\n"; my $command = "$warptool -setskyfiletoupdate -warp_id $warp_id -skycell_id $skycell->{skycell_id}"; $command .= " -set_label $rlabel" if $rlabel; @@ -527,4 +544,5 @@ } } elsif ($chips_ready and $skycell->{data_state} eq 'update' and $skycell->{state} ne 'update') { + print "Chips are ready for warp $warp_id $skycell->{skycell_id} but run state is $skycell->{state}. Changing state to update.\n"; my $command = "$warptool -updaterun -warp_id $warp_id -set_state update"; $command .= " -set_label $rlabel" if $rlabel; @@ -566,9 +584,8 @@ my $skycell_id = $skycell->{skycell_id}; - # if ($metadata->{state} eq 'goto_cleaned' and $metadata->{data_state} eq 'full') { if ($metadata->{state} eq 'goto_cleaned') { # cleanup must not be running. Set state to update. If the skycell is not 'full' it will be - # set to be updated the next time this script is invoked - print "changing state of diffRun $diff_id from goto_cleaned to update. skycell $skycell_id is full\n"; + # set to be updated the next time this script is invoked. + print "Changing state of diffRun $diff_id from goto_cleaned to update.\n"; my $command = "$difftool -updaterun -set_state update -diff_id $diff_id"; $command .= " -set_label $rlabel" if $rlabel; @@ -589,5 +606,5 @@ # check the state of the template stack my $command = "$stacktool -sumskyfile -stack_id $skycell->{stack2}"; - my $stack = runToolAndParseExpectOne($command, $verbose); + my $stack = runToolAndParseExpectOne($command, $parse_verbose); my_die("failed to find stackSumSkyfile for stack_id $skycell->{stack2}", $PS_EXIT_UNKNOWN_ERROR) if !$stack; @@ -597,8 +614,9 @@ return $PSTAMP_GONE; } + print " Stack $skycell->{stack2} is ready.\n"; # now check the warp $command = "$warptool -warped -warp_id $skycell->{warp1} -skycell_id $skycell_id"; - my $warp = runToolAndParseExpectOne($command, $verbose); + my $warp = runToolAndParseExpectOne($command, $parse_verbose); my_die("failed to find warpSkyfile for warpRun $skycell->{warp1} skycell_id $skycell_id", $PS_EXIT_UNKNOWN_ERROR) if !$warp; @@ -615,7 +633,8 @@ } # warps are ready fall through and queue the diff update + print " Warp $skycell->{warp1} $skycell_id is ready.\n"; } elsif ($diff_mode eq $IPP_DIFF_MODE_WARP_WARP) { my $command = "$warptool -warped -warp_id $skycell->{warp1} -skycell_id $skycell_id"; - my $warp1 = runToolAndParseExpectOne($command, $verbose); + my $warp1 = runToolAndParseExpectOne($command, $parse_verbose); my_die("failed to find warpSkyfile for warpRun $skycell->{warp1} skycell_id $skycell_id", $PS_EXIT_UNKNOWN_ERROR) if !$warp1; @@ -637,5 +656,5 @@ } $command = "$warptool -warped -warp_id $skycell->{warp2} -skycell_id $skycell_id"; - my $warp2 = runToolAndParseExpectOne($command, $verbose); + my $warp2 = runToolAndParseExpectOne($command, $parse_verbose); my_die("failed to find warpSkyfile for warpRun $skycell->{warp2} skycell_id $skycell_id", $PS_EXIT_UNKNOWN_ERROR) if !$warp2; @@ -661,5 +680,5 @@ # check the state of the input stack my $command = "$stacktool -sumskyfile -stack_id $skycell->{stack1}"; - my $stack1 = runToolAndParseExpectOne($command, $verbose); + my $stack1 = runToolAndParseExpectOne($command, $parse_verbose); my_die("failed to find stackSumSkyfile for stack_id $skycell->{stack1}", $PS_EXIT_UNKNOWN_ERROR) if !$stack1; @@ -671,5 +690,5 @@ # check the state of the template stack $command = "$stacktool -sumskyfile -stack_id $skycell->{stack2}"; - my $stack2 = runToolAndParseExpectOne($command, $verbose); + my $stack2 = runToolAndParseExpectOne($command, $parse_verbose); my_die("failed to find stackSumSkyfile for stack_id $skycell->{stack2}", $PS_EXIT_UNKNOWN_ERROR) if !$stack2; @@ -684,5 +703,5 @@ # check the state of the input stack my $command = "$stacktool -sumskyfile -stack_id $skycell->{stack1}"; - my $stack = runToolAndParseExpectOne($command, $verbose); + my $stack = runToolAndParseExpectOne($command, $parse_verbose); my_die("failed to find stackSumSkyfile for stack_id $skycell->{stack1}", $PS_EXIT_UNKNOWN_ERROR) if !$stack; @@ -695,5 +714,5 @@ # now check the template warp $command = "$warptool -warped -warp_id $skycell->{warp2} -skycell_id $skycell_id"; - my $warp = runToolAndParseExpectOne($command, $verbose); + my $warp = runToolAndParseExpectOne($command, $parse_verbose); my_die("failed to find warpSkyfile for warpRun $skycell->{warp2} skycell_id $skycell_id", $PS_EXIT_UNKNOWN_ERROR) if !$warp; @@ -736,4 +755,6 @@ print "skipping $command\n"; } + } else { + print " diff is in update state\n"; } @@ -763,5 +784,7 @@ my $now = DateTime->now->mjd; my $dtime_tool = (DateTime->now->mjd - $start_tool) * 86400.; - print "Time to run $program: $dtime_tool\n"; + if ($dtime_tool > 0.1) { + print "Time to run $program: $dtime_tool\n"; + } my $buf = join "", @$stdout_buf; @@ -777,5 +800,7 @@ my $dtime_parse = (DateTime->now->mjd - $start_parse) * 86400.; - print "Time to parse results from $program: $dtime_parse\n"; + if ($dtime_parse > 0.1) { + print "Time to parse results from $program: $dtime_parse\n"; + } return $results; @@ -848,5 +873,5 @@ if (!$dsRun_state) { my $command = "$chiptool -listrun -chip_id $stage_id"; - my $data = runToolAndParse($command, $verbose); + my $data = runToolAndParse($command, $parse_verbose); my $chipRun = $data->[0]; $dsRun_state = $chipRun->{dsRun_state}; @@ -862,5 +887,5 @@ foreach my $c (@$components) { my $command = "$magicdstool -destreakedfile -magic_ds_id $magic_ds_id -component $c"; - my $dsfile = runToolAndParseExpectOne($command, $verbose); + my $dsfile = runToolAndParseExpectOne($command, $parse_verbose); if (!$dsfile) { my_die("failed to find magicDSFile for ${stage}Run $stage_id $c", $PS_EXIT_UNKNOWN_ERROR); Index: branches/eam_branches/ipp-20140904/tools/heather/check_stsci0.pl =================================================================== --- branches/eam_branches/ipp-20140904/tools/heather/check_stsci0.pl (revision 37804) +++ branches/eam_branches/ipp-20140904/tools/heather/check_stsci0.pl (revision 37804) @@ -0,0 +1,6 @@ +@num = ('00','01','02','03','04','05','06','07','08','09','10','11','12','13','14','15','16','17','18','19'); + +foreach $numb (@num) { + $thing = `ssh stsci$numb 'df -h | grep stsci$numb.0'`; + print $thing; +} Index: branches/eam_branches/ipp-20140904/tools/heather/check_stsci1.pl =================================================================== --- branches/eam_branches/ipp-20140904/tools/heather/check_stsci1.pl (revision 37804) +++ branches/eam_branches/ipp-20140904/tools/heather/check_stsci1.pl (revision 37804) @@ -0,0 +1,6 @@ +@num = ('00','01','02','03','04','05','06','07','08','09','10','11','12','13','14','15','16','17','18','19'); + +foreach $numb (@num) { + $thing = `ssh stsci$numb 'df -h | grep stsci$numb.1'`; + print $thing; +} Index: branches/eam_branches/ipp-20140904/tools/heather/check_stsci2.pl =================================================================== --- branches/eam_branches/ipp-20140904/tools/heather/check_stsci2.pl (revision 37804) +++ branches/eam_branches/ipp-20140904/tools/heather/check_stsci2.pl (revision 37804) @@ -0,0 +1,6 @@ +@num = ('00','01','02','03','04','05','06','07','08','09','10','11','12','13','14','15','16','17','18','19'); + +foreach $numb (@num) { + $thing = `ssh stsci$numb 'df -h | grep stsci$numb.2'`; + print $thing; +} Index: branches/eam_branches/ipp-20140904/tools/heather/queuesliceslap2.pl =================================================================== --- branches/eam_branches/ipp-20140904/tools/heather/queuesliceslap2.pl (revision 37803) +++ branches/eam_branches/ipp-20140904/tools/heather/queuesliceslap2.pl (revision 37804) @@ -4,9 +4,11 @@ -#$label = "LAP.ThreePi.20130717"; -$label = "LAP.ThreePi.20130717.pole3"; -#$label = "LAP.ThreePi.20130717.pole3.CNP"; - +$label1 = "LAP.ThreePi.20130717"; +$label2 = "LAP.ThreePi.20130717.pole3"; +$label3 = "LAP.ThreePi.20130717.pole3.CNP"; +#heather done bad things to segments and slices... +@labels = ($label1, $label2, $label3); @segments = (0,1,2,3,4,5,6,7,8); +#@segments=(0,1); @segments1 = ( -35, -22,-10,2,14,26,38,50,62); @segments2 = ( -22,-10,2,14,26,38,50,62,75); @@ -14,5 +16,8 @@ @slices = ("00","02","04","06","08","10","12","14","16","18","20","22"); - +#@slices = ("00","02","04"); + +foreach $label (@labels) { + foreach $slice (@slices) { @@ -21,6 +26,6 @@ foreach $num (@segments) { - $dvodb = "LAP.PV2.20130717"; - $minidvodb_group = "LAP.PV2.slice$slice.$segments[$num].cam"; + $dvodb = "LAP.PV2.V2.20130717"; + $minidvodb_group = "LAP.PV2.V2.slice$slice.$segments[$num].cam"; $dvodb_label = "$minidvodb_group"; $ra_min = $slice*15.0; @@ -37,5 +42,6 @@ } -$com2 = "addtool -definebyquery -label $label -set_dvodb $dvodb -set_minidvodb_group LAP.PV2.pole.cam -set_minidvodb -set_label LAP.PV2.pole.cam -stage cam -uncensored -ra_min 0 -ra_max 360 -decl_min 75 -decl_max 90 -dbname gpc1"; +$com2 = "addtool -definebyquery -label $label -set_dvodb $dvodb -set_minidvodb_group LAP.PV2.V2.pole.cam -set_minidvodb -set_label LAP.PV2.V2.pole.cam -stage cam -uncensored -ra_min 0 -ra_max 360 -decl_min 75 -decl_max 90 -dbname gpc1"; print "$com2\n"; `$com2`; +} Index: branches/eam_branches/ipp-20140904/tools/heather/queuesliceslap2skycal.pl =================================================================== --- branches/eam_branches/ipp-20140904/tools/heather/queuesliceslap2skycal.pl (revision 37803) +++ branches/eam_branches/ipp-20140904/tools/heather/queuesliceslap2skycal.pl (revision 37804) @@ -4,8 +4,8 @@ -#$label = "LAP.ThreePi.20130717"; -#$label = "LAP.ThreePi.20130717.pole3"; -$label = "LAP.ThreePi.20130717.pole3.CNP"; - +$label1 = "LAP.ThreePi.20130717"; +$label2 = "LAP.ThreePi.20130717.pole3"; +$label3 = "LAP.ThreePi.20130717.pole3.CNP"; +@labels = ($label1,$label2,$label3); @segments = (0,1,2,3,4,5,6,7,8); @segments1 = ( -35, -22,-10,2,14,26,38,50,62); @@ -14,5 +14,5 @@ @slices = ("00","02","04","06","08","10","12","14","16","18","20","22"); - +foreach $label (@labels) { foreach $slice (@slices) { @@ -21,6 +21,6 @@ foreach $num (@segments) { - $dvodb = "LAP.PV2.20130717"; - $minidvodb_group = "LAP.PV2.slice$slice.$segments[$num].skycal"; + $dvodb = "LAP.PV2.V2.20130717"; + $minidvodb_group = "LAP.PV2.V2.slice$slice.$segments[$num].skycal"; $dvodb_label = "$minidvodb_group"; $ra_min = $slice*15.0; @@ -36,6 +36,8 @@ } + } -$com2 = "addtool -definebyquery -label $label -set_dvodb $dvodb -set_minidvodb_group LAP.PV2.pole.skycal -set_minidvodb -set_label LAP.PV2.pole.skycal -stage skycal -uncensored -ra_min 0 -ra_max 360 -dec_min 75 -dec_max 90 -dbname gpc1"; +$com2 = "addtool -definebyquery -label $label -set_dvodb $dvodb -set_minidvodb_group LAP.PV2.V2.pole.skycal -set_minidvodb -set_label LAP.PV2.V2.pole.skycal -stage skycal -uncensored -ra_min 0 -ra_max 360 -dec_min 75 -dec_max 90 -dbname gpc1"; print "$com2\n"; `$com2`; +} Index: branches/eam_branches/ipp-20140904/tools/iasc/getit =================================================================== --- branches/eam_branches/ipp-20140904/tools/iasc/getit (revision 37803) +++ branches/eam_branches/ipp-20140904/tools/iasc/getit (revision 37804) @@ -45,5 +45,5 @@ chdir $final_destdir or die "failed to chdir to $final_destdir"; -system 'perl /data/ippc18.0/home/ipp/iasc.v2/rename.stamp.outputs.pl *'; +system 'perl /data/ippc18.0/home/ipp/iasc/rename.stamp.outputs.pl *'; Index: branches/eam_branches/ipp-20140904/tools/stask_chip.py =================================================================== --- branches/eam_branches/ipp-20140904/tools/stask_chip.py (revision 37804) +++ branches/eam_branches/ipp-20140904/tools/stask_chip.py (revision 37804) @@ -0,0 +1,112 @@ +#!/usr/bin/env python + +import os +import time +import sys +import subprocess + +def emit_output(p): + (sout, serr) = p.communicate(None) + if len(sout): + print >> sys.stdout, sout, + sys.stdout.flush() + +def wait(n): + global n_completed + while (len(procs) >= n): + t = time.time() + for p in procs: + if p.poll() is not None: + print >> sys.stderr, '%.2f -- %s is done, status %d' % (t-times[p], cmds[p].rstrip(), p.returncode) + emit_output(p) + n_completed=n_completed+1 + busy[prochost[p]] = 0 + del prochost[p] + del cmds[p] + del times[p] + procs.remove(p) + elif t-times[p] > timeout: + print >> sys.stderr, '%s timed out!' % cmds[p].rstrip() + p.terminate() + emit_output(p) + busy[prochost[p]] = 0 + del prochost[p] + del cmds[p] + del times[p] + procs.remove(p) + time.sleep(1.0) + time.sleep(0.05) + +def free_host(): + for host in hosts: + if busy[host] != 1: + return host + return None + +task_file = sys.argv[1] +hosts_file = sys.argv[2] +nprocs = int(sys.argv[3]) # Oversubscription parameter: number of processors to devote to a task + +hosts = [] +tasks = [] +nname = -1 # prefix of host specification +procs = [] # process, defined as task on a host +cmds = {} # command assocated with process +prochost = {} # host associate with process +times = {} +busy = {} +timeout = 2*3600 +n_active = 0 # number of active jobs we expect to have running should equal Nhosts * Ncore_per_host / nprocs + + +n_completed=0 +start_time = time.time() + +with open(task_file) as f: + for task in f: + if task[0] != '#': + tasks.append(task.rstrip()) + +with open(hosts_file) as f: + for host in f: + if host[:2] == 'mu': + cores = int(12 / nprocs) # Use nprocs to reduce the core subscription + nname = 6 + elif host[:2] == 'ml': + cores = 16 + nname = 5 + elif host[:2] == 'mp': + cores = 8 + nname = 5 + else: + cores = 1 + nname = -1 + for core in range(cores): + h = "%s-%d" % (host[:nname], core) + hosts.append(h) + busy[h] = 0 + n_active = n_active + 1 + +print 'running %d commands at a time from %d tasks in %s' \ + % (n_active, len(tasks), task_file) + +for task in tasks: + h = free_host() + busy[h] = 1 + command = '/usr/bin/ssh -n -q %s \"cd /scratch3/watersc1/; %s\"' % (h[:nname], task) + print 'executing ', command + sys.stdout.flush() + try: + p = subprocess.Popen(command, bufsize=-1, shell=True, stdin=None, stdout=subprocess.PIPE, stderr=subprocess.STDOUT) + except: + print >> sys.stderr, "Exception! %s" % command + procs.append(p) + prochost[p] = h + cmds[p] = command + times[p] = time.time() + wait(n_active) + +wait(1) + +t = time.time() +print >> sys.stderr, 'all done. %d tasks completed in %2.f seconds (%.2f per second)' % (n_completed, t-start_time, n_completed/(t-start_time))