Subs::UvotImages (version 0.0)


package Subs::UvotImages;
##############################################################################
#
# DESCRIPTION: Extract transform and plot images
#
# HISTORY:
# HISTORY: $Log: UvotImages.pm,v $
# HISTORY: Revision 1.44  2005/02/08 18:21:32  apsop
# HISTORY: Add new types for uvot filter specific images and exp maps.  Replace Filename::corresponding method with file name hack for now.  Fix bug which was preventing the production of raw images.
# HISTORY:
# HISTORY: Revision 1.43  2005/01/21 04:26:43  apsop
# HISTORY: Change to use file names stored in cat file rather than contructing them.
# HISTORY:
# HISTORY: Revision 1.42  2005/01/12 17:29:23  apsop
# HISTORY: Change weightfile param in uvotdetect to expfile param.
# HISTORY:
# HISTORY: Revision 1.41  2005/01/08 04:11:00  apsop
# HISTORY: Change test for existence of event list to be more forgiving.
# HISTORY:
# HISTORY: Revision 1.40  2005/01/07 20:34:21  apsop
# HISTORY: No longer necessary to use the CENTER method of swiftxform or a small
# HISTORY: value of ncell for uvotmodmap.
# HISTORY:
# HISTORY: Revision 1.39  2004/12/03 13:38:11  apsop
# HISTORY: Turn back on calling uvotmodmap.
# HISTORY:
# HISTORY: Revision 1.38  2004/10/13 01:38:31  apsop
# HISTORY: Disable uvotmodmap for the time being, and fix bug in getting plot file names.
# HISTORY:
# HISTORY: Revision 1.37  2004/10/12 16:26:36  apsop
# HISTORY: Filename changes and a few other tweaks to make module run in complete pipeline.
# HISTORY:
# HISTORY: Revision 1.36  2004/09/15 22:26:24  apsop
# HISTORY: Remove the binning factor from uvot event list file names
# HISTORY:
# HISTORY: Revision 1.35  2004/09/01 14:42:21  apsop
# HISTORY: Disable level 2/3 processing for now.
# HISTORY:
# HISTORY: Revision 1.34  2004/08/31 20:51:02  apsop
# HISTORY: Implemented UVOT level 2 and 3 products.
# HISTORY:
# HISTORY: Revision 1.33  2004/08/18 17:27:58  apsop
# HISTORY: Yet another fix for the timing keywords.
# HISTORY:
# HISTORY: Revision 1.32  2004/07/19 16:03:24  apsop
# HISTORY: Add version number entry to comment field.
# HISTORY:
# HISTORY: Revision 1.31  2004/07/12 13:46:09  apsop
# HISTORY: Fix problems with timing keywords
# HISTORY:
# HISTORY: Revision 1.30  2004/06/15 16:01:04  apsop
# HISTORY: Better trapping of problem uvot exposures.
# HISTORY:
# HISTORY: Revision 1.29  2004/06/02 16:11:19  apsop
# HISTORY: Change name of catalogue extension.
# HISTORY:
# HISTORY: Revision 1.28  2004/05/04 16:33:01  dah
# HISTORY: Bug fixes for case where no initial image file exists.
# HISTORY:
# HISTORY: Revision 1.27  2004/04/30 16:16:40  dah
# HISTORY: Use new FITSfile::keywords() method to speed up image extraction
# HISTORY:
# HISTORY: Revision 1.26  2004/04/16 20:21:18  dah
# HISTORY: Begin using embedded history records
# HISTORY:
# HISTORY: Revision 1.25  2004/04/16 19:53:18  dah
# HISTORY: Begin using embedded history records
# HISTORY:
#
# VERSION: 0.0
# 
##############################################################################


use Subs::Sub;
use Subs::Images;
use Util::Xanadu;
use Subs::UvotNames;

@ISA = ("Subs::Images");
use strict;

sub new {
    my $proto=shift;
    my $self=$proto->SUPER::new();

    $self->{DESCRIPTION}="Extracting, merging and plotting images for UVOT";

    return $self;
}

##################
# METHODS:
##################

sub body {
    my $self=shift;

    #####################################
    # extract images from the event data
    #####################################
    $self->extract_images();

    $self->process_raw();

    $self->image_raw_to_det();

    $self->image_raw_to_sky('uvot');

    $self->create_exposure_maps();

    # create summed images/exposure maps
    $self->sum_images();

    # detect sources and determine magnitudes
    $self->detect_sources();

    $self->plot_images();


} # end of body method


#############################################################################
# Extract images from filtered event data
#############################################################################
sub extract_images {
    my $self=shift;

    my $filterCodes = $Subs::UvotNames::filterCodes;
    my $filterNames = $Subs::UvotNames::filterParams;

    my @exp_keywords = ('TSTART', 'TSTOP', 'TELAPSE', 'ONTIME', 'LIVETIME', 'EXPOSURE', 'DEADC');
    my @ds_keywords = ('DSUNI1', 'DSVAL1', 'DSTYP1', 'DSREF2', 'DSUNI2', 'DSVAL2', 'DSTYP2');

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();

    $log->entry("Extracting images from uvot event data as required");

    my $catFile = $filename->get('hk', 'uvot', 'ct', '*');

    unless( -f $catFile ){
      $log->entry("No catalogue file, so cannot extract images.");
      return;
    }

    my $catFits = Util::FITSfile->new($catFile, 'EXPCATALOG');
    unless ($catFits->keyword('EXTNAME')) {
      $log->error(2, "File $catFile does not have an EXPCATALOG extension.");
      return;
    }

    my @filterId = $catFits->cols('FILTER')->table();
    my @gtiName = $catFits->cols('GTINAME')->table();
    my @expRef = $catFits->cols('EXPREF')->table();
    my @efile = $catFits->cols('EFILEREF')->table();
    my @ifile = $catFits->cols('IFILEREF')->table();

    my @binning = $catFits->cols('BINLVL')->table();
    my @hasImage = $catFits->cols('Image')->table();
    my @hasEvents = $catFits->cols('Event')->table();
    my @startTime = $catFits->cols('ESTART')->table();
    my @stopTime = $catFits->cols('ESTOP')->table();

    my @x0 = $catFits->cols('EV_X0')->table();
    my @y0 = $catFits->cols('EV_Y0')->table();
    my @xsize = $catFits->cols('EV_XSIZ')->table();
    my @ysize = $catFits->cols('EV_YSIZ')->table();

    my @imx0 = $catFits->cols('IM_X0')->table();
    my @imy0 = $catFits->cols('IM_Y0')->table();
    my @imxsize = $catFits->cols('IM_XSIZ')->table();
    my @imysize = $catFits->cols('IM_YSIZ')->table();

    my @winx0 = $catFits->cols('WINDOWX0')->table();
    my @winy0 = $catFits->cols('WINDOWY0')->table();
    my @windx = $catFits->cols('WINDOWDX')->table();
    my @windy = $catFits->cols('WINDOWDY')->table();
    my @binx = $catFits->cols('BINX')->table();
    my @biny = $catFits->cols('BINY')->table();

    my $catRows = $catFits->nrows();

    ###########################
    # fire up ftools
    ###########################
    my $fextract = Util::Ftool->new('fextract');
    my $tabedit = Util::HEAdas->new('ftedit');

    my %needsSorting;
    my $tmpFile = 'Uvot_raw_images.tmp';

    for(my $exp=0; $exp<$catRows; $exp++){
      if(!$binning[$exp] || $binning[$exp] eq 'INDEF'){
	$log->entry("Exposure row $exp has an undefined or zero binning value." .
		    "  Skipping extraction of this exposure.");
	next;
      }
      next if $hasEvents[$exp] eq 'b0';
      ##################################################################
      # Only make images for those exposures that don't already have one
      # or if the event window is not completely contained in the image
      # window
      ##################################################################
      unless($hasImage[$exp] eq 'b0'){
	my $extract =  ($x0[$exp]-$xsize[$exp]/2) < ($imx0[$exp]-$imxsize[$exp]/2) 
	            || ($y0[$exp]-$ysize[$exp]/2) < ($imy0[$exp]-$imysize[$exp]/2)
		    || ($x0[$exp]+$xsize[$exp]/2) > ($imx0[$exp]+$imxsize[$exp]/2)
		    || ($y0[$exp]+$ysize[$exp]/2) > ($imy0[$exp]+$imysize[$exp]/2);
	next unless $extract;
      }

      my $expref = $expRef[$exp];
      $binning[$exp] =~ s/\s//g;

      my $evtFile = $efile[$exp];
      my $mode = $$filterCodes[$filterId[$exp]] . 'po';
      my $imageFile = $filename->get('rawimage', 'uvot', $$filterCodes[$filterId[$exp]], 0);
      my $extname = $$filterCodes[$filterId[$exp]] . sprintf('%08d', int($startTime[$exp])) . 'E';
				    
      my $evtFits = Util::FITSfile->new($evtFile, 'WINDOW');
      my @windExpRef = $evtFits->cols('EXPREF')->table();
       
      ###################################
      # Update catalogue and event files
      ###################################
      $tabedit->params({infile => $catFile.'[1]',
			row => $exp+1});

      $tabedit->params({column => 'IFILEREF',
			value => $imageFile})
              ->run();

      $tabedit->params({column => 'IMENAME',
			value => $extname})
              ->run();

      $tabedit->params({column => 'Image',
			value => 128})
              ->run();

      for(my $irow=0; $irow<@windExpRef; $irow++){
	if($windExpRef[$irow]==$expref){
	  $tabedit->params({infile => $evtFile.'[WINDOW]',
			    row => $irow+1,
			    column => 'IMAGEXT',
			    value => 1})
                  ->run();

	}
      }

      ##############################
      # Set up to extract the image
      ##############################
      $evtFits->ext('EVENTS');
      $evtFits->specs('[gtifilter(\'[' . $gtiName[$exp] . ']\')]' .
		      '[bin RAWX=' . $winx0[$exp] . ':' .
		                     ($winx0[$exp]+$windx[$exp]-1) . ':' .
		                     $binx[$exp] . ',' .
		           'RAWY=' . $winy0[$exp] . ':' .
		                     ($winy0[$exp]+$windy[$exp]-1) . ':' .
                                     $biny[$exp] . ']');

      $log->entry("Extracting raw image from $evtFile, gti " . $gtiName[$exp] .' '. `date`);

      if( -f $imageFile ){
	$needsSorting{$imageFile} = 1;
      }else{
	Util::Ftool->new('fimgcreate')
	           ->params({bitpix => 8,
			      naxes => '0',
			     datafile => 'none',
			     outfile => $imageFile})
                   ->run();
      }
      $evtFits->append_to($imageFile);

      my $imageFits = Util::FITSfile->new($imageFile, 0);
      $imageFits->ext($imageFits->nhdus()-1);

      ###############################################################
      # Fix up the keywords.  You'd think this would be easy ...
      ###############################################################
      my %keywords = $imageFits->keywords();
      my @keynames = keys %keywords;

      $imageFits->begin_many_keywords();

      $imageFits->keyword('EXTNAME', $extname);

      my ($cdelt1, $cdelt2) = ($imageFits->keyword('CDELT1'), $imageFits->keyword('CDELT2'));

      $imageFits->keyword('WINDOWX0', 0.5 + (0.5-$imageFits->keyword('CRPIX1'))*$cdelt1 + 
			  $imageFits->keyword('CRVAL1') );
      $imageFits->keyword('WINDOWY0', 0.5 + (0.5-$imageFits->keyword('CRPIX2'))*$cdelt1 + 
			  $imageFits->keyword('CRVAL2') );
      $imageFits->keyword('WINDOWDX', $imageFits->keyword('NAXIS1')*$cdelt1 );
      $imageFits->keyword('WINDOWDY', $imageFits->keyword('NAXIS2')*$cdelt2 );
      $imageFits->keyword('BINX', $cdelt1);
      $imageFits->keyword('BINY', $cdelt2);

      $imageFits->keyword('CMPUTHRS', -1000);
      $imageFits->keyword('CMPOTHRS', 1000);
      $imageFits->keyword('CMPCNTMN', 'F');
      $imageFits->keyword('BLANK', -1);

      $imageFits->keyword('TSTART', $startTime[$exp]);
      $imageFits->keyword('TSTOP', $stopTime[$exp]);
      $imageFits->keyword('TELAPSE', $stopTime[$exp]-$startTime[$exp] );

      my $first = Util::Date->new($startTime[$exp]);
      my $last  = Util::Date->new($stopTime[$exp]);

      $imageFits->keyword('DATE-OBS', $first->date().'T'.$first->time() );
      $imageFits->keyword('DATE-END', $last->date().'T'.$last->time() );

      foreach my $exp_key (@exp_keywords){
	for(my $index=1; $index<=$catRows; $index++){
	  my $oldname = substr($exp_key, 0, 6) . $index;
	  $imageFits->keyword( "-$oldname", ' ' ) if grep /${oldname}/, @keynames;

	  if($index==$expref){
	    $imageFits->keyword( $exp_key, $evtFits->keyword($oldname) );
	  }
	}
      }

      foreach my $ds_key (@ds_keywords){
	for(my $index=1; $index<=$catRows; $index++){
	  my $oldname = $index==1 ? '' : $index;
	  $oldname .= $ds_key;
	  $imageFits->keyword( "-$oldname", ' ' ) if grep /${oldname}/, @keynames;
	}
      }

      $imageFits->keyword( 'HDUCLAS1', 'IMAGE' );
      $imageFits->keyword( 'HDUCLAS2', 'TOTAL' );

      $imageFits->end_many_keywords();

      my $livetime = $imageFits->keyword('ONTIME');
      my $deadtime = $imageFits->keyword('DEADC');
      $livetime = $livetime*$deadtime if $deadtime;

      $imageFits->keyword( 'LIVETIME', $livetime);
      $imageFits->keyword( 'EXPOSURE', $livetime);
    }

    #############################################################
    # Sort the extensions to be in the right order, if necessary
    #############################################################

    my $fappend = Util::Ftool->new("fappend");
    foreach my $sortFile (keys %needsSorting) {
      unlink $tmpFile;
      $fextract->params({'infile'  => $sortFile . '[0]', 'outfile' => $tmpFile})
	->run();

      my @hduNames;
      my $sortFits = Util::FITSfile->new($sortFile);
      my $nhdus = $sortFits->nhdus();
      for(my $hdu=1; $hdu<$nhdus; $hdu++){
	$sortFits->ext($hdu);
	push @hduNames, $sortFits->keyword('EXTNAME');
      }
      my @sortedHduNames = sort @hduNames;
      next if join('%',@sortedHduNames) eq join('%',@hduNames);

      foreach my $extName (@sortedHduNames) {
	next unless $extName;
	$fappend->params({'infile'  => $sortFile . "[$extName]",
			  'outfile' => $tmpFile})
	  ->run();
      }
      rename $tmpFile, $sortFile;
    }
    unlink $tmpFile;

} # end extract images method

###############################################################################
# uvot specific processing of the raw images
###############################################################################
sub process_raw {
    my $self=shift;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();

    $log->entry("Processing raw uvot images");

    my $badpix = $filename->fetch_cal('badpix', 'u');
    my $ubadpix = Util::HEAdas->new('uvotbadpix');
    $ubadpix->params({'infile' => $badpix});

    my $umodmap = Util::HEAdas->new('uvotmodmap');

    my $flat = $filename->fetch_cal('flat', 'u');
    my $uflatfield = Util::HEAdas->new('uvotflatfield');
    $uflatfield->params({flatfile => $flat});

    $log->entry("setting {RA,DEC,PA}_PNT in all raw images");
    my $ra = $jobpar->read("ra");
    my $dec = $jobpar->read("dec");
    my $roll = $jobpar->read("roll");

    ##############################
    # iterate over the UVOT images
    ##############################
    foreach my $rawFile ($filename->get('rawimage', 'uvot', '*', '*')) {
      {
 	my $imageFits = Util::FITSfile->new($rawFile);
	for my $i (1 .. $imageFits->nhdus - 1) {
	  $imageFits->ext($i);
	  $imageFits->keyword('RA_PNT', $ra);
	  $imageFits->keyword('DEC_PNT', $dec);
	  $imageFits->keyword('PA_PNT', $roll);
 	}
      }

      my $badFile = $filename->corresponding('rawimage', 'badimage', $rawFile);
      my $corrFile = $filename->corresponding('rawimage', 'corrimage', $rawFile);
      my $mod8File = 'mod.tmp';

      $log->entry("running uvotbadpix on $rawFile");
      $ubadpix->params({infile => $rawFile,
                        badpixlist => $badpix,
			outfile => $badFile,
			chatter => 3})
	      ->run();

      $log->entry("running uvotmodmap on $rawFile");
      $umodmap->params({infile => $rawFile,
                        badpixfile => $badFile,
			outfile => $mod8File,
			nsig => 3,
			ncell => 16, # TODO: keeps run times for 2048^2 images to minutes
			            # instead of hours
			chatter => 3})
              ->run();

      # uflatfield will go here
      $uflatfield->params({
            infile => $mod8File,
	    infile => $rawFile,
            outfile => $corrFile})
            ->run();

      unlink($mod8File);
    }
}

###############################################################################
# convert raw coordinate images to detector coordinates
###############################################################################
sub image_raw_to_det {
    my $self=shift;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();

    $log->entry("Converting raw coordinate images to detector coordinates");

    #########################################
    # create the swiftxform tool
    #########################################
    my $swiftxform = Util::HEAdas->new('swiftxform');

    my $teldef = $filename->fetch_cal('teldef', 'uvot');
    $swiftxform->params({
			 teldeffile => $teldef,
			 to         => 'DET',
			 attfile    => 'NONE',
			 aberration => 'no',
			 seed       => $procpar->read('seed'),
			 chatter    => 4,
			})
                ->is_script(1);

    ######################
    # loop over _GRISM_ files
    #   since only they are converted to DET coordinates
    ######################
    foreach my $rawFile ($filename->get('corrimage', 'uvot', 'g*', '*') ) {

      #########################################
      # determine the DET coordinate filename
      #########################################
      my $detFile = $filename->corresponding("corrimage", "detimage", $rawFile);
      $log->entry("converting $rawFile to $detFile");

      ######################
      # do the conversion
      ######################
      $swiftxform->params({
            infile   => $rawFile,
            outfile  => $detFile,
            })
         ->run();

    } # end of loop over files

} # end of image_raw_to_det method


###############################################################################
# convert raw coordinate images to sky coordinates
###############################################################################
sub image_raw_to_sky {
    my $self=shift;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();

    $log->entry("Converting raw coordinate images to sky coordinates");


    ######################################################
    # make sure there is an attitude file
    ######################################################
    my $attitude = $filename->get('attitude', 's');
    unless(-f $attitude) {
        $log->entry("No attitude data available - can't make sky images");
        return;
    }

    #########################################
    # create the swiftxform tool
    #########################################
    my $swiftxform = Util::HEAdas->new('swiftxform')->is_script(1);

    my $teldef = $filename->fetch_cal('teldef', 'uvot');
    $swiftxform->params({
             teldeffile => $teldef,
             to         => 'SKY',
             attfile    => $attitude,
             ra         => $jobpar->read('ra'),
             dec        => $jobpar->read('dec'),
             method     => 'AREA',
             aberration => 'no',
             seed       => $procpar->read('seed'),
             chatter    => 4,
             });

    ######################
    # loop over files
    ######################
    foreach my $rawFile ($filename->get('corrimage', 'uvot', '*', '*') ) {

      #########################################
      # determine the DET coordinate filename
      #########################################
      my $skyFile = $rawFile;
      $skyFile =~ s/_cr\./_sk\./;
      $log->entry("converting $rawFile to $skyFile");

      ######################
      # do the conversion
      ######################
      $swiftxform->params({
            infile   => $rawFile,
            outfile  => $skyFile,
            })
         ->run();

    } # end of loop over files

} # end of image_raw_to_sky method


###############################################################################
# generate exposure maps
###############################################################################
sub create_exposure_maps {
    my $self=shift;

    my $log     = $self->log();
    my $filename= $self->filename();
    my $procpar = $self->procpar();
    my $jobpar  = $self->jobpar();

    $log->entry("Creating exposure maps");


    ######################################################
    # make sure there is an attitude file
    ######################################################
    my $attitude = $filename->get('attitude', 's');
    unless(-f $attitude) {
        $log->entry("No attitude data available - can't make exposure maps");
        return;
    }

    #########################################
    # create the swiftxform tool
    #########################################
    my $uexpmap = Util::HEAdas->new('uvotexpmap');

    my $teldef = $filename->fetch_cal('teldef', 'uvot');
    $uexpmap->params({
             attfile    => $attitude,
             teldeffile => $teldef,
             method     => 'MEANFOV', # TODO: is SHIFT/ADD in effect?
             attdelta   => 100, # TODO: make a parameter
             aberration => 'no',
             chatter    => 4,
             });

    ######################
    # loop over files
    ######################
    foreach my $skyFile ($filename->get('filterimg', 'uvot', '*', '*') ) {

      # exclude grism files
      next if $filename->is_grism($skyFile, 'skyimage');

      #########################################
      # determine the exposure map filename
      #########################################
      my $expFile = $skyFile;
      $expFile =~ s/_sk\./_ex\./;
      $log->entry("building exposure map $expFile for $skyFile");

      my $badFile = $skyFile;
      $badFile =~ s/_sk\./_bp\./;

      ######################
      # do the conversion
      ######################
      $uexpmap->params({
            infile => $skyFile,
            badpixfile => $badFile,
            outfile => $expFile,
            })
         ->is_script(1)
         ->run();

    } # end of loop over files

} # end of create_exposure_maps method



sub sum_images {

    my $self=shift;

    my $log     = $self->log();
    my $filename= $self->filename();
    my $procpar = $self->procpar();
    my $jobpar  = $self->jobpar();

    $log->entry("Creating summed images");


    ################################
    # Create tool for summing images
    ################################
    my $usumexp = Util::HEAdas->new('uvotimsum');

    my $tmpFile = 'sum.tmp';

    ########################################
    # clean up any pre-existing summed files
    ########################################
    foreach my $type (qw(skyimage expimage)) {

       foreach my $imageFile ($filename->get($type, 'uvot', '', '*')) {

          $log->entry("removing obsolete summed image $imageFile");
          unlink($imageFile);
       }
    }

    ######################
    # loop over files
    ######################
    my %filter_to_sum = ('filterimg' => 'skyimage',
			 'filterexp' => 'expimage');

    foreach my $type (keys %filter_to_sum) {

       foreach my $imageFile ($filename->get($type, 'uvot', '*', '*')) {

          my (undef, $filter_code) = $filename->parse($imageFile, $type);
	  $log->entry("Filter code for $imageFile is $filter_code.");

          next if not $filter_code;
      
          #########################################
          # determine the summed image filename
          # by removing the filter code from sw\d{9}uff(sk|ex).img
          #########################################
          my $sumFile = $filename->get($filter_to_sum{$type}, 'uvot', '', 0);
          $log->entry("summing $imageFile to $sumFile");

          my $filter = $filename->filter_name_for_code($filter_code);

          ######################
          # do the summation
          ######################
          $usumexp->params({
                infile => $imageFile,
                outfile => $tmpFile,
                })
             ->is_script(1)
             ->run();

          my $keys = "[col "
                 . " #EXTNAME='$filter';"
                 . " #FILTER='$filter';"
                 . ']';

          if (not -e $sumFile) {
             Util::HEAdas->new('ftcopy')
                ->params({
                      infile => $tmpFile . "[col #EXTNAME='$filter']",
                      outfile => $sumFile,
                      copyall => 'yes',
                   })
                ->run();
          }
          else {
             # append this result to the sum file
             Util::HEAdas->new('ftappend')
                ->params({
                      infile => "$tmpFile+1[col #EXTNAME='$filter']",
                      outfile => $sumFile,
                   })
                ->run();
          }

           unlink($tmpFile);
        } # end of loop over files

    } # end of loop over types

} # end of process_sky method



sub detect_sources {

    my $self=shift;

    my $log     = $self->log();
    my $filename= $self->filename();
    my $procpar = $self->procpar();
    my $jobpar  = $self->jobpar();

    $log->entry("Detecting sources");


    ########################################
    # clean up any pre-existing source lists
    ########################################
    foreach my $catFile ($filename->get('srclist', 'uvot', '', '*')) {
       unlink($catFile);
    }
      
    ###############################################################
    # Create tools for detecting sources and determining magnitudes
    ###############################################################
    my $udetect = Util::HEAdas->new('uvotdetect');
    my $umag = Util::HEAdas->new('uvotmag');
    my $color = $filename->fetch_cal('color', 'u');
    my $cntcor = $filename->fetch_cal('cntcor', 'u');
    $umag->params({
        zerofile => $color,
        coinfile => $cntcor,
        filter => 'default',
        });

    my $tmpFile = 'sources.tmp';
    my $catFile = undef;

    #############################
    # loop over summed sky images
    #############################
    foreach my $skyFile ($filename->get('skyimage', 'uvot', '', '*')) {

        $log->entry("detect: skyFile '$skyFile'");

        my $expFile = $filename->corresponding('skyimage', 'expimage', $skyFile);

        if (not $expFile or not -f $expFile) {
            $log->entry("No exposure data for $skyFile - assuming constant exposure");
            $expFile = 'NONE';
        }

        if (not $catFile) {
            $catFile = $filename->corresponding('skyimage', 'srclist', $skyFile);
        }

        #########################################
        # run uvotdetect on each extension
        #########################################
        my $fitsFile = Util::FITSfile->new($skyFile, 0);
        my $nhdus = $fitsFile->nhdus();

        for (my $hdu0 = 1; $hdu0 < $nhdus; ++$hdu0) {

            $fitsFile->ext($hdu0);
            my $filter = $fitsFile->keyword('EXTNAME');

            my $weightFile = $expFile;
            if ($expFile ne 'NONE') {
               $weightFile .= "[$filter]";
            }

            if (-e $tmpFile) {
               unlink($tmpFile);
            }

            $udetect->params({
                 infile => "$skyFile+$hdu0",
                 outfile => $tmpFile,
                 expfile => $weightFile,
                 threshold => 2,  # TODO: make parameter
                 })
              ->is_script(1)
              ->run();

            # run uvotmag
            if (-e $tmpFile) {
               $umag->params({
                    infile => $tmpFile,
                    })
                 ->is_script(1)
                 ->run();
            }

            if (not -e $catFile) {
               Util::HEAdas->new('ftcopy')
                  ->params({
                        infile => $tmpFile . "[col #EXTNAME='$filter']",
                        outfile => $catFile,
                        copyall => 'yes',
                     })
                  ->run();
            }
            else {
               # append this result to the sum file
               Util::HEAdas->new('ftappend')
                  ->params({
                        infile => "$tmpFile+1[col #EXTNAME='$filter']",
                        outfile => $catFile,
                     })
                  ->run();
            }
        } # end of loop over extensions

    } # end of loop over files

    unlink($tmpFile) if -e $tmpFile;

} # end of detect_sources method


################################################################################
#
################################################################################
sub plot_images {
    my $self = shift;

    my $log     =$self->log();
    my $filename=$self->filename();
    my $procpar =$self->procpar();
    my $jobpar  =$self->jobpar();

    my $max_image_dimen=1024;

    $log->entry("Plotting all images");
    my $ximage = Util::Xanadu->new("ximage");

    my $type;
    foreach $type ("detimage", "skyimage") {

        my $plot_type = $type;
        $plot_type =~ s/image/plot/;

        ######################################
        # loop over summed images of this type
        ######################################
        my $file;
        foreach my $file ($filename->get($type, 'uvot', '', '*')) {

            my $fitsFile = Util::FITSfile->new($file, 0);
            my $nhdus = $fitsFile->nhdus();

            my $grid='';
            if($type =~ /^sky/ ) { $grid = 'grid' }

            ##############################################
            # for each extension (filter) create a graphic
            ##############################################
            for(my $hdu=1; $hdu<$nhdus; $hdu++) {

               my @commands;
               push @commands, ('cpd /vgif', 'cey 2000');

               $fitsFile->ext($hdu);
               my $filter = $fitsFile->keyword('EXTNAME');
               my $filter_code = $filename->filter_code_for_name($filter);

               my $plotFile = $filename->get($plot_type, 'uvot', $filter_code, 0);

               ########################################
               # determine if we need to bin the image
               ########################################
               my $dimen = $fitsFile->keyword('NAXIS1');
               my $bin=1;
               for($bin=1; $max_image_dimen*$bin<$dimen; $bin++) {}

               $log->entry("Plotting ${file}\[$hdu\] in $plotFile binned by $bin");

               push @commands, ("read/fits/rebin=$bin \{${file}\[$hdu\]\}",
                                'smooth',
                                'disp',
                                $grid,
                                'scale');

               push @commands, 'exit';
               $ximage->script(@commands)->run();
               rename 'pgplot.gif', $plotFile;
            }
        
        }
    }

} # end of plot_images method

1;