Subs::BATImages (version $)

package Subs::BATImages;
##############################################################################
#
# DESCRIPTION: Accumulate detector plane images in the various combinations
#		that are needed.
#
#		Port of Craig Markwardt's BAT::dpi, image modules to SDC.
#
#
# HISTORY: $Log: BATImages.pm,v $
# HISTORY: Revision 1.9  2005/04/06 15:41:05  apsop
# HISTORY: Change to using CALDB for cal parameters.
# HISTORY:
# HISTORY: Revision 1.8  2005/03/07 20:24:01  apsop
# HISTORY: Change detimage type to dpimage.
# HISTORY:
# HISTORY: Revision 1.7  2005/02/10 00:08:35  apsop
# HISTORY: Added message identifiers.
# HISTORY:
# HISTORY: Revision 1.6  2005/02/08 19:15:46  apsop
# HISTORY: Fix problem with BAT image names which was causing images not to be produced.
# HISTORY:
# HISTORY: Revision 1.5  2004/12/31 01:34:09  apsop
# HISTORY: Comment out 4 channel images, as they are not official products.  Change 1chan images to proper file name.
# HISTORY:
# HISTORY: Revision 1.4  2004/12/06 00:58:17  apsop
# HISTORY: Fix typo - info() should be entry()
# HISTORY:
# HISTORY: Revision 1.3  2004/12/05 23:41:23  apsop
# HISTORY: Change pcodeimg error message to an info message
# HISTORY:
# HISTORY: Revision 1.2  2004/11/16 14:24:26  apsop
# HISTORY: Only try to process files that exist.
# HISTORY:
# HISTORY: Revision 1.1  2004/11/08 18:40:39  apsop
# HISTORY: Module for creating BAT image products.
# HISTORY:
# HISTORY: Revision 1.1  2004/09/24 20:51:19  wiegand
# HISTORY: Initial revision
# HISTORY:
#
# VERSION: $Revision: 1.9 $
#
#
##############################################################################


use Subs::Sub;
use Util::Parfile;
use Util::Ftool;
use Util::HEAdas;
use Util::BATCave;
use Util::CoreTags;
use Util::SwiftTags;


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

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

	$self->{DESCRIPTION}= 'Processing BAT images';

	return $self;
}

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

sub body {

	my $self = shift;

	$self->detector_plane_images;

	$self->sky_images;

	$self->partial_coding_maps;

	$self->postage_stamp_images;

	$self->glom_images;

}


sub detector_plane_images
{
	my ($self) = @_;

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

	$log->entry('Accumulating BAT detector plane images');


	my @config = (
		{ key => 'preburst',
			events => 'evshps',
			gti => 'GTI_BKG1',
			mode => 'evpb',
		},
		{ key => 'preslew',
			events => 'evshps',
			gti => 'GTI_TOT',
			mode => 'evps',
		},
		{ key => 'postslew',
			events => 'evshas',
			gti => 'NONE',
			mode => 'evas',
		},
	);


	my %info = (
		qmap   => $filename->get('qualcal', 'bat', 'cb'),
	);


	foreach my $c (@config) {

		my @events = $filename->getExisting(
				'unfiltered', 'bat', $c->{events});
		next if not @events;

		$info{files} = \@events;
		$info{gti} = Util::BATCave::get_gti($self, $c->{gti});

		$info{energybins} = Util::BATCave::chan1;
		$info{outfile} = $filename->get('dpimage', 'bat',
				$c->{mode}, 0),
		$self->dpiMake(\%info);

#		$info{energybins} = Util::BATCave::chan4;
#		$info{outfile} = $filename->get('dpimage', 'bat',
#				$c->{mode} . '4chan', 0),
#		$self->dpiMake(\%info);
	}


} # end of body method



# BAT::dpi->make

sub dpiMake
{
	my ($self, $info) = @_;

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

	my $dpifile = $info->{outfile};
	my $erange = $info->{energybins} || '-';
	my $qmap = $info->{qmap} || 'NONE';
	my $gti = $info->{gti} || 'NONE';

	my $inlist = join(',', @{ $info->{files} });

	unlink($dpifile);

	my $batbinevt = Util::HEAdas->new('batbinevt')
			->params({
				infile     => $inlist,
				outfile    => $dpifile,
				outtype    => 'DPI',
				timedel    => 0,
				timebinalg => 'uniform',
				energybins => $erange,
				detmask    => $qmap,
				ecol       => 'ENERGY',
				weighted   => 'NO',
				outunits   => 'COUNTS',
				gtifile    => $gti,
				clobber    => 'yes',
				chatter    => 3,
				})
			->run;

	return if $batbinevt->had_error;

	# Note: this case is not fatal, *IF* there were no overlapping
	# good times between the input file and the good time interval
	# file.  In that case batbinevt does not create an output file.

	if (not -f $dpifile) {
        $log->entry("warning: batbinevt did not create DPI $dpifile");
    }

}



sub sky_images
{
	my ($self) = @_;

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

	$log->entry('Computing sky images ...');

	my %common = (
		teldef => 'CALDB',
		attfile => $filename->get('attitude', 's'),

		corrections => 'autocollim,flatfield,ndets,pcode,maskwt',
		aperture => 'CALDB',
		qmap => $filename->get('qualcal', 'bat', 'cb'),

		bat_z => $jobpar->read('bat_z'),
		origin_z => $jobpar->read('bat_origin_z'),

		# Partial coding threshold.  Minimum detector plane exposure
		# threshold.
		pcodethresh => 0.05, # Fraction of 1.0
	);

	# TODO: check for %common files or return


	my @config = (
		{ interval => 'preburst',
			mode => 'evpb',
		},
		{ interval => 'preslew',
			mode => 'evps',
		},
		{ interval => 'postslew',
			mode => 'evas',
		},
	);

	foreach my $c (@config) {

		### my @chan_ranges = qw(1chan 4chan);

		# can truncate to 1chan if $run_for_speed

		my $done = 0;
		my @files;

		### foreach my $chan (@chan_ranges) {

			# TODO: file types?
			my $mode = $c->{mode};

			my $dpifile = $filename->get('dpimage', 'bat', $mode);
			next if not -f $dpifile;

			my $imgfile = $filename->get('skyimage', 'bat', $mode, 0);

			my $bkgfile = 'NONE';
			if ($c->{interval} eq 'preslew') {
				# Pre-slew interval uses preburst background
				$bkgfile = $filename->get('dpimage', 'bat', 'evpb');
			}

			my %info = (
				%common,
				dpifile => $dpifile,
				imgfile => $imgfile,
				bkgfile => $bkgfile,
			);

			$self->imageSky(\%info);
			$done = 1;
		### }

		if ($done) {
			$log->entry("... $c->{interval} done");
		}
	}

}


# BAT::image->sky

sub imageSky
{
	my ($self, $info) = @_;

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

	my $imgfile = $info->{imgfile};

	my $bkgfile = $info->{bkgfile} || 'NONE';
	my $pcodemap = $info->{pcodemap} || 'NO';

	unlink($imgfile);

	# XXX NOTE: clobber=NO is a workaround for a bug in the build 9
	# batfftimage task.  Since we unlink($imgfile) just above, this is
	# equivalent to clobber=yes

	my $batbinevt = Util::HEAdas->new('batfftimage')
			->params({
				infile     => $info->{dpifile},
				outfile    => $imgfile,
				attitude   => $info->{attfile},
				bkgfile    => $bkgfile,
				bat_z      => $info->{bat_z},
				origin_z   => $info->{origin_z},
				teldef     => $info->{teldef},
				pcodethresh=> $info->{pcodethresh},
				corrections=> $info->{corrections},
				detmask    => $info->{qmap},
				aperture   => $info->{aperture},
				pcodemap   => $pcodemap,
				clobber    => 'no',
				chatter    => 3,
				})
			->run;

}


sub partial_coding_maps
{
	my ($self) = @_;

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

	$log->entry('Computing partial coding maps...');

	my %common = (
		teldef => 'CALDB',
		attfile => $filename->get('attitude', 's'),

		# since pcodemap == imageSky(bkgfile=NONE, pcodemap=YES)
		bkgfile   => 'NONE',
		pcodemap  => 'YES',

		corrections => 'autocollim,flatfield,ndets,pcode,maskwt',
		aperture => 'CALDB',
		qmap => $filename->get('qualcal', 'bat', 'cb'),

		bat_z => $jobpar->read('bat_z'),
		origin_z => $jobpar->read('bat_origin_z'),

		# Partial coding threshold.  Minimum detector plane exposure
		# threshold.
		pcodethresh => 0.05, # Fraction of 1.0
	);


	# Only the pre-slew and post-slew segments are meaningful for the
	# partial coding map (during the slew doesn't make sense)

	my @config = (
		{ interval => 'preslew',
			mode => 'evps',
		},
		{ interval => 'postslew',
			mode => 'evas',
		},
	);

	foreach my $c (@config) {

		my $done = 0;
		# Only one partial coding map per interval
		foreach my $chan (qw(1chan)) {

			my $mode = $c->{mode} . $chan;

			# TODO: file types?
			my $dpifile   = $filename->get('dpimage', 'bat', $mode);
			next if not -f $dpifile;
			my $pcodefile = $filename->get('expimage', 'bat', $c->{mode}, 0);

			my %info = (
				%common,
				dpifile => $dpifile,
				imgfile => $pcodefile,
			);

			$self->imageSky(\%info);
			$done = 1;
		}   

		if ($done) {
			$log->entry("... $c->{interval} done");
		}
	}

	my $pcodeimg = $filename->get('expimage', 'bat', 'evps');
	if (not -f $pcodeimg) {
		$log->entry('no BAT preslew pcodeimg, unable to set pcode');
		return;
	}

	my $ra = $jobpar->read('burst_ra');
	my $dec = $jobpar->read('burst_dec');
	my $value = $self->pcodeFromImage($pcodeimg, $ra, $dec);
	if (defined($value)) {
		if (undef) {
			$jobpar->set({
					bat_pcode => $value,
					});
		}
		$log->entry("determined BAT pcode => $value");
	}
	else {
		$log->error(1, 'unable to set BAT pcode');
	}
}


sub sky2pix
{
	my ($image, $ra, $dec) = @_;

	my $sky2xy = Util::Ftool->new('sky2xy')
			->params({
				infile => $image,
				xsky => $ra,
				ysky => $dec,
				})
			->run;

	return if $sky2xy->had_error;

	# grab xpix, ypix
	my $skypars = Util::Parfile->new('sky2xy.par');
	my $xpix = $skypars->read('xpix');
	my $ypix = $skypars->read('ypix');

	return [ $xpix, $ypix ];
}



# BAT::pcode->fromimage / BAT::imageutils
sub pcodeFromImage
{
	my ($self, $path, $ra, $dec) = @_;

	my $pixpos = sky2pix($path, $ra, $dec);
	return if not $pixpos;

	my ($xpix, $ypix) = map { sprintf('%d', $_) } @$pixpos;

	my $fimgpar = Util::Ftool->new('fimgpar')
			->params({
				fitsfile => $path,
				pixel => "$xpix,$ypix",
				})
			->run;
	return if $fimgpar->had_error;

	my $imgpars = Util::Parfile->new('fimgpar.par');
	my $undef = $imgpars->read('undef');
	return if $undef =~ m/^y/i;

	my $value = $imgpars->read('outvalue');
	return $value;
}



sub postage_stamp_images
{
	my ($self) = @_;

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

	$log->entry('Extracting postage stamp images...');

	my $ra = $jobpar->read('burst_ra');
	my $dec = $jobpar->read('burst_dec');

	my $postsize = 100;

	my @config = (
		{ interval => 'preburst',
			mode => 'evpb',
		},
		{ interval => 'preslew',
			mode => 'evps',
		},
		{ interval => 'postslew',
			mode => 'evas',
		},
	);


	foreach my $c (@config) {

		# Only one partial coding map per interval
		foreach my $chan (qw(1chan)) {

			my $mode = $c->{mode} . $chan;

			my $imgfile = $filename->get('skyimage', 'bat', $mode);
			next if not -f $imgfile;

			my $postfile = $filename->get('postimg', 'bat', $mode, 0);
			unlink($postfile);

			my $pixpos = sky2pix($imgfile, $ra, $dec);
			if (not $pixpos) {
				$log->error([ 1, BAT_TASK_ERROR ],
					"unable to make postage stamp image for $imgfile");
				next;
			}

			my ($xpix, $ypix) = @$pixpos;

			my $xstart = int($xpix - $postsize / 2);
			my $xstop = $xstart + $postsize;

			my $ystart = int($ypix - $postsize / 2);
			my $ystop = $ystart + $postsize;

			if ($xstart < 1) { $xstart = 1 };
			if ($ystart < 1) { $ystart = 1 };

			my $copier = Util::HEAdas->new('ftcopy')
				->params({
						infile  => $imgfile . "[$xstart:$xstop,$ystart:$ystop]",
						outfile => $postfile,
					})
				->run;
		}
	}
}



sub glomImages
{
	my ($self, $outfile, @f) = @_;

	my $log = $self->log;

	foreach my $f (@f) {
		if (-e $outfile) {
			my $appender = Util::HEAdas->new('ftappend')
					->params({
							infile => $f . '[0]', # TODO: which extension(s)?
							outfile => $outfile,
						})
					->run;
			if ($appender->had_error) {
				$log->error([ 2, HEATOOL_ERROR ],
					"unable to append $f to $outfile");
			}
		}
		else {
			rename($f, $outfile)
				or $log->error([ 2, BAD_OUTPUT ],
					"unable to rename $f => $outfile [$!]");
		}
	}
}


sub glom_images
{
	my $self = shift;

	my $log = $self->log;

	$log->entry('this is where images would be packaged up');
}




1;