Changeset 34041 for branches/meh_branches/ppstack_test/tools/dlapstacks

Timestamp:

Jun 19, 2012, 5:24:19 PM (14 years ago)

Author:

mhuber

Message:

merging latest r34040 trunk changes to branch

Location:

branches/meh_branches/ppstack_test

Files:

• . (modified) (1 prop)
• tools/dlapstacks (modified) (11 diffs)

branches/meh_branches/ppstack_test/tools/dlapstacks

@@ -4 +4 @@
 use warnings;
 
-# given an RA and DEC lookup find the LAP stack images
-
-# XXX: consider making these options these could be made options I guess
-
+# given an RA and DEC or a skycell_id find the LAP stack images containing them
+# and display using ds9
+
+# XXX: consider accepting a string of parameters to be passed to ds9
+# add --fit set ds9 to fit to frame
+# XXX if label is supplied by the user and it is like MD% 
+# guess at the tess_id
+# allow mask and variance images to be displayed as well
+
+my $label   = 'LAP.ThreePi.20110809';
 my $tess_id = 'RINGS.V3';
-my $label = 'LAP.ThreePi.20110809';
 
 use Getopt::Long qw( GetOptions :config auto_help auto_version gnu_getopt );
@@ -15 +20 @@
 use DBI;
 use PS::IPP::Config 1.01 qw( :standard );
-
-my $ipprc =  PS::IPP::Config->new();
+use File::Temp qw( tempfile tempdir );
+use File::Basename;
+use IPC::Cmd 0.36 qw( can_run run );
+
+my $ipprc = PS::IPP::Config->new();
 my $dbh = getDBHandle();
 
-
 my ($requested_filter, $no_display, $skycell_id, $verbose, $ra, $dec, $convolved);
-
+my ($singleframe, $fitToFrame, $deasin, $save_temps, $data_group);
+
+# XXX: get this from the ipprc
+my $temproot = '/local/ipp/tmp';
+my $outdir;
 
 GetOptions(
-    'filter|f=s'    => \$requested_filter,
-    'convolved'     => \$convolved,
-    'no-display|n'  => \$no_display,
-    'skycell_id=s'  => \$skycell_id,
-    'verbose|v'     => \$verbose,
+    'filter|f=s'    => \$requested_filter,  # restrict to a single filter 
+    'skycell_id|s=s'=> \$skycell_id,        # specificly sset skycell_id of interest
+    'convolved|c'   => \$convolved,         # display convolved stack images instead of unconvolved
+    'fit'           => \$fitToFrame,        # scale the images to fit the ds9 frame
+    'single'        => \$singleframe,       # set ds9 to singleframe mode
+    'deasin|d'      => \$deasin,            # run ppstamp to remove the ASIN scaling from the stacks
+    'label|l=s'     => \$label,             # use a different label (careful!)
+    'data_group=s'  => \$data_group,        # select stacks by data_group
+    'tess_id|t=s'   => \$tess_id,           # set tess_id
+    'outdir|o=s'    => \$outdir,            # if deasin save de-asined images to given directory
+    'no-display|n'  => \$no_display,        # don't run ds9 just print the nebulous pathnames
+    'save-temps'    => \$save_temps,        # save temporary files
+    'verbose|v'     => \$verbose,           # print lots of information
 ) or pod2usage( 2 );
 
@@ -36 +55 @@
 $ra = shift;
 $dec = shift;
+
+pod2usage( -msg => 'if ra is suplied <dec> must be supplied as well', -exitval => 2)
+    if (defined $ra and ! defined $dec);
+
+if ($skycell_id and !($skycell_id =~ /^skycell/)) {
+    if (length($skycell_id) < 3) {
+        # two character skycell_id NM should probably be 0NM
+        $skycell_id = "0$skycell_id";
+    }
+    $skycell_id = "skycell.$skycell_id";
+}
+
+my $missing_tools;
+my $ppstamp;
+my $ds9;
+if ($deasin) {
+    $ppstamp = can_run('ppstamp') or (warn "Can't find ppstamp") and $missing_tools = 1;
+}
+if (!$no_display) {
+    $ds9 = can_run('ds9') or (warn "Can't find ds9") and $missing_tools = 1;
+}
+if ($missing_tools) {
+    # warn("Can't find required tools.");
+    exit($PS_EXIT_CONFIG_ERROR);
+}
 
 my $extension;
@@ -44 +88 @@
 }
 
-# run whichimage to find a list of skycells containing the given coordinates
-my $whichimage_output = `whichimage --listchipcoords --tess_id $tess_id $ra $dec`;
-
-print $whichimage_output if $verbose;
-
-my @results = split "\n", $whichimage_output;
-
-die "whichimage returned no output\n" if scalar @results == 0;
-
-my $query = "SELECT skycell_id, filter, stack_id, quality, path_base
+# don't allow MD label unless it's a refstack.
+# Nightly stacks have too many images
+if (($data_group && ($data_group =~ /MD/)) or ($label =~ /^MD/)) {
+    unless ($data_group or uc($label) =~ /REF/) {
+        die "cannot use this program for label $label without data_group because that will most likely select nightly stacks and we will match an unreasonably large number of matches. Talk to Bill\n";
+    }
+    # if tess_id hasn't been set to an MD field, guess that it's a V3 stack
+    unless ($tess_id =~ /MD/) {
+        my $tag = $data_group ? $data_group : $label;
+        my $field = substr $tag, 0, 4;
+        $tess_id = "$field.V3";
+        print "Setting tess_id to $tess_id\n";
+    }
+}
+
+my $whichimage_output;
+if (!$skycell_id or (defined $ra and defined $dec)) {
+    # run whichimage to find a list of skycells containing the given coordinates and the 
+    # image coordinates of the requested sky coordinates
+    $whichimage_output = `whichimage --listchipcoords --tess_id $tess_id $ra $dec`;
+
+    print "Output from whichimage:\n$whichimage_output\n" if $verbose;
+
+    my @results = split "\n", $whichimage_output;
+
+    die "whichimage returned no output\n" if scalar @results == 0;
+} else {
+    # we were supplied a skycell create a fake whichimage_output that can be parsed below
+    $whichimage_output = "0 0 dummy $skycell_id -1 -1\n";
+}
+
+my $query = "SELECT skycell_id, filter, stack_id, quality, state, path_base
     FROM stackRun join stackSumSkyfile USING(stack_id) 
-    WHERE stackRun.state ='full' AND label = '$label' AND skycell_id = ?";
+    WHERE stackRun.state ='full' AND skycell_id = ?";
+if ($data_group) {
+    $query .= " AND data_group = '$data_group'";
+} else {
+    $query .= " AND label = '$label'";
+}
 
 if ($requested_filter) {
@@ -69 +140 @@
 
 # @files will contain the list of skycells found in order of skycell_id, filter.
-# These will be displayed in this order so that the image frames are in order of increasing bandpass
+# The images will be displayed in this order so that the image frames are in order of increasing bandpass
 my @files = ();
 
@@ -81 +152 @@
 foreach (split "\n", $whichimage_output) {
     chomp;
+
     my ($radeg, $decdeg, undef, $skycell, $x, $y) = split " ";
-    print "$skycell $radeg $decdeg\n" if $verbose;
+
+    print "skycell: $skycell RA (deg): $radeg DEC (deg): $decdeg (X,Y): ($x, $y)\n" if $verbose;
+
     $stmt->execute($skycell) or die $stmt->errstr;
     my $files_found = 0;
     while (my $f = $stmt->fetchrow_hashref()) {
         next if $f->{quality};
+        next if $f->{state} ne 'full';
+
         my $filter = $f->{filter};
         my $stack_id = $f->{stack_id};
         my $path_base = $f->{path_base};
+
         # save the skycell coordinates so we can pan efficiently below
         $f->{X} = $x;
         $f->{Y} = $y;
 
-        # Find where this filter goes in the array
+        # Find where this filter goes in the array of files
         my $i = $index{$filter};
         if (!defined $i) {
@@ -100 +177 @@
             next;
         }
-        # if multiple skycells match the coordinates bump the index
+        # if multiple skycells match the coordinates bump the index by $num_filters
         $i += $skycells * $num_filters;
         $files[$i] = $f;
@@ -112 +189 @@
 
 if (!$total_files) {
-    print "no files found\n";
+    print "no matching stacks found\n";
     exit 1;
 }
@@ -118 +195 @@
 
 if (!$no_display) {
-    # estimate of image size for the case where we have 5 frames
-    my $imagesize = 680;
-    # my $panargs = "-pan to $ra $dec wcs fk5";
+    my $imagesize;
+    if ($singleframe) {
+        $imagesize = 1024;
+    } else {
+        # estimate of image size for the case where we have 5 frames
+        $imagesize = 680;
+    }
+    my $tempdir;
+    if (!$outdir) {
+        $tempdir = tempdir ("$temproot/dlapstacks.temp.XXXX", CLEANUP => !$save_temps);
+    }
+    if ($deasin) {
+        if ($outdir) {
+            print "Saving deasined files to $outdir\n" if $verbose;
+            # user specified output directory save the files there
+            if (!-d $outdir) {
+                mkdir $outdir or die "failed to mkdir $outdir\n";
+            }
+        } else {
+            # otherwise use a temporary directory
+            $outdir = $tempdir;
+            print "deasined files will be saved to $outdir\n" if $save_temps;
+        }
+    }
+    # Build the ds9 argument list
     my @args;
+    my $panargs = "";
     foreach my $f (@files) {
         next if ! $f;
-        my $filename = $f->{path_base} . ".$extension";
-        my $resolved = $ipprc->file_resolve($filename);
-        next if !$resolved;
-        my $panargs = "-pan to $f->{X} $f->{Y}";
-
+        my $pathname = $f->{path_base} . ".$extension";
+        my $filename = $ipprc->file_resolve($pathname);
+        next if !$filename;
+
+        # prefix file name with filter
+        # this will either be the ppstamp output filename or we
+        # will put a symlink to it in tempdir
+        my $base = basename($pathname);
+        my $filt = substr $f->{filter}, 0, 1;
+        $base = "$filt.$base";
+        if ($deasin) {
+            # run ppstamp to undo the asin scaling. The output from ppstamp is linear F32 images
+
+            # drop the .fits for the output path for ppstamp
+            my $destbase = $base;
+            $destbase =~ s/\.fits//g;
+
+
+            my $command = "$ppstamp -wholefile -file $filename $outdir/$destbase";
+            my ( $success, $error_code, $full_buf, $stdout_buf, $stderr_buf ) =
+                run(command => $command, verbose => $verbose);
+            unless ($success) {
+                print STDERR join "", @$stderr_buf;
+                $error_code = (($error_code >> 8) or $PS_EXIT_PROG_ERROR);
+                die("Unable to perform ppstamp: $error_code\n");
+            }
+
+            # new filename for ds9
+            $filename = "$outdir/$base";
+        } else {
+            my $nebname = $filename;
+            $filename = "$tempdir/$base";
+            symlink $nebname, $filename or die "failed to symlink $nebname";
+        }
+
+        if ($f->{X} > 0) {
+            $panargs = "-pan to $f->{X} $f->{Y}";
+        } else {
+            $panargs = "";
+        }
+
+        # The following is an attempt put the filter name in the image using a region with text a label.
         my $xlabel = $f->{X} + $imagesize / 2;
         my $ylabel = $f->{Y} + $imagesize / 2;
-        # The following is an attempt put the filter name in the image using a region with text a label.
+
         # I couldn't make this work. The circle would appear but without the text
         # my $label = " -regions command  'circle($xlabel,$ylabel,1)' # text={ $f->{filter} }";
@@ -140 +277 @@
         # $panargs .= $label;
 
-        push @args, ($resolved, $panargs);
-    }
-    my $ds9cmd = "ds9 -geom 1024x1024 @args";
+        # I used to pan each frame but it turns out to be quicker to wait until the last image is
+        # loaded. pan thatone and then matching the wcs for the frames.
+        # Also that allows V0 and V3 skycells to have east to the left orientation
+        # push @args, ($filename, $panargs);
+        push @args, $filename;
+    }
+    my $ds9cmd = "$ds9 -geom 1024x1024 @args $panargs";
+    $ds9cmd .= " -single" if $singleframe;
+    $ds9cmd .= " -zoom to fit" if $fitToFrame;
+    $ds9cmd .= " -match frames wcs";
+    $ds9cmd .= " -frame first" if $singleframe;
     print "$ds9cmd\n" if $verbose;
     system $ds9cmd;
 } else {
+    # no-display option Just print out some information
+    print "\n" if $verbose;
+    print "Skycell ID        stack_id filter   path_base\n";
+    print "------------------------------------------------------------------------------------------\n";
     foreach my $f (@files) {
         next if !$f;
-        print "$f->{skycell_id} $f->{stack_id} $f->{filter} $f->{path_base}\n";
+        print "$f->{skycell_id}  $f->{stack_id}   $f->{filter}  $f->{path_base}\n";
     }
 }
@@ -173 +322 @@
 }
 
+__END__
+
+=pod
+=head1 NAME
+dlapstacks - display stack images from the PS1 first Grand Reprocessing (LAP)
+
+=head1 SYNOPSIS
+    
+    dlapstacks <ra> <decl> [options]
+    dlapstacks --skycell_id <skycell_id> [<ra> <decl>] [options]
+
+
+=head1 DESCRIPTION
+
+Query the gpc1 database to find stacks containing a given RA and DEC or for a specific skycell_id.
+Display the images using ds9. If coordinates are supplied the display is centered there.
+
+=head1 OPTIONS
+
+=over 4
+
+=item * -s, --skycell_id <skycell_id>
+
+The skycell_id to display. Note if skycell_id does not begin with "skycell", "skycell" is prepended to the supplied value before querying the database.
+
+
+=item * <ra>   
+
+Right Ascension. Format: decimal degrees or HH:MM:SS
+Required unless skycell_id is supplied
+
+
+=item * <decl>
+
+Declination. Format: decimal degrees or DD:MM:SS 
+Required unless skycell_id is supplied
+
+=item * -f, --filter <filter>
+
+Only display images for the given filter. The SQL wild card character '%' is appended to single character values.
+
+=item * -c, --convolved
+
+Display convolved stacks. Default is to display unconvolved stack images.
+
+=item * --fit
+
+Instruct ds9 to fit the image to the frame.
+
+=item * --single
+
+Set ds9 to single frame mode
+
+=item * -d, --deasin
+
+Display images tha have had the ASIN scaling removed.
+
+=item *-o, --outdir <outdir>
+
+If --deasin is selected the unscaled images are preserved in the provided directory.
+The directory created if it does not exist.
+
+=item * -l, --label <label>
+
+Chose stacks with the provided label. The default is LAP.ThreePi.20110809.
+NOTE: If the label begins with MD and --data_group is not provided the label must contain "ref" or "REF". This prevents
+users from selecting the hundreds of nightly stack images. 
+
+=item * --data_group <data_group>
+
+Select stacks with the provided data_group. By default stacks are selected by label.
+
+
+=item * -t, --tess_id <tess_id>
+
+Set the tessellation that skycells are selected from. Default is RINGS.V3
+unless a label or data group that begins with MDnm is supplied. In that case the default is
+MDnm.V3
+
+=item * -n, --no-display
+
+Look up the stack images using the provided parameters but rather than displaying them with ds9 print out some information about the stacks.
+
+=item * --save-temps
+
+If --deasin is supplied and --outdir is not supplied do not delete temporary images. (For debugging).
+
+=item * -v, --verbose
+
+Print information about the progress of the program.
+
+=cut