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.27 2014/02/27 09:22:47 apsop
# HISTORY: Don't use uat attitude file, at request of BAT Team/C.Markwardt.
# HISTORY:
# HISTORY: Revision 1.26 2012/08/30 07:04:48 apsop
# HISTORY: sub plot_images: Note source of trigtime. (Info not currently used,
# HISTORY: done for symmetry with eg. UvotProduct::AddTrigtime.)
# HISTORY:
# HISTORY: Revision 1.25 2012/01/12 06:52:03 apsop
# HISTORY: Changes going to proc3.15.03
# HISTORY:
# HISTORY: 2011-09-19 JRG as apsop: When deleting HDUNAME kwd, use deleteall
# HISTORY: which doesn't complain if it's missing as delete does.
# HISTORY:
# HISTORY: Revision 1.24 2011/01/20 17:22:02 apsop
# HISTORY: Added code to use UVOT attitude file if available
# HISTORY: Modified plotting subroutine to create plots consistent with format
# HISTORY: used in UVOT plots
# HISTORY: Added routine rearrange_images so that extesnion BAT_IMAGE_TOT is the
# HISTORY: first one
# HISTORY:
# HISTORY: Revision 1.23 2008/05/15 20:20:31 apsop
# HISTORY: Run batimgstatpos on imag status telemetry.
# HISTORY:
# HISTORY: Revision 1.22 2007/10/23 20:52:17 apsop
# HISTORY: Use burst pos to make postage stamp images, instead of pointing pos.
# HISTORY:
# HISTORY: Revision 1.21 2007/01/31 21:15:10 apsop
# HISTORY: Do not make postage stamp images outside the FOV. Improve and enable postage stamp plots.
# HISTORY:
# HISTORY: Revision 1.20 2006/10/29 19:30:19 apsop
# HISTORY: Use BAT_DPI_TOT ext name instead of BAT_IMAGE_TOT for detector images;
# HISTORY: First version of plotting code (not run).
# HISTORY:
# HISTORY: Revision 1.19 2006/09/10 20:07:22 apsop
# HISTORY: Add new distfile parameter for build 19.
# HISTORY:
# HISTORY: Revision 1.18 2006/08/02 19:47:30 apsop
# HISTORY: Change failure of sky2xyto a level 1 error.
# HISTORY:
# HISTORY: Revision 1.17 2006/08/01 20:57:35 apsop
# HISTORY: Trap the case where batbinevt does not produce an output image due to an
# HISTORY: empty filtered file. Continue processing.
# HISTORY:
# HISTORY: Revision 1.16 2006/07/06 18:16:09 apsop
# HISTORY: Change ext name of 1chan image to BAT_IMAGE_TOT. Test for presence of
# HISTORY: BAT_CATALOG ext in 1chan image before copying.
# HISTORY:
# HISTORY: Revision 1.15 2006/05/10 18:42:37 apsop
# HISTORY: Better trapping of empty source list, and keep ps images from extracting
# HISTORY: beyond the edges of the input image.
# HISTORY:
# HISTORY: Revision 1.14 2006/05/10 15:26:55 apsop
# HISTORY: Trap case where batcelldetect finds no sources.
# HISTORY:
# HISTORY: Revision 1.13 2005/11/08 17:21:51 apsop
# HISTORY: Change calls for qualcal and skyimage filenames. Append BAT_CATALOG
# HISTORY: extension to preslew image.
# HISTORY:
# HISTORY: Revision 1.12 2005/07/06 20:36:55 apsop
# HISTORY: Implemented BAT position refinement.
# HISTORY:
# HISTORY: Revision 1.11 2005/05/16 14:15:33 apsop
# HISTORY: Handle missing input files more gracefully when creating images. Pass
# HISTORY: CALDB for ebounds parameter when creating mask tagged light curves.
# HISTORY:
# HISTORY: Revision 1.10 2005/05/13 19:38:02 apsop
# HISTORY: Implemented product agreements between HEASARC and BAT team.
# HISTORY:
# 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
# HISTORY: 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.27 $
#
#
##############################################################################
use Subs::Sub;
use Util::Parfile;
use Util::Ftool;
use Util::HEAdas;
use Util::BATCave;
use Util::CoreTags;
use Util::SwiftTags;
use Util::SimpleFITS;
use Astro::FITS::CFITSIO qw(:constants);
use Subs::UvotProduct;
@ISA = ("Subs::Sub");
use strict;
sub new {
my $proto=shift;
my $self=$proto->SUPER::new();
$self->{DESCRIPTION}= 'Processing BAT images';
$self->{POSCAT} = 'posrefine.cat';
$self->{POSREG} = 'posrefine.reg';
return $self;
}
##################
# METHODS:
##################
sub body {
my $self = shift;
$self->detector_plane_images;
$self->sky_images;
# $self->partial_coding_maps;
$self->postage_stamp_images;
$self->refine_position;
$self->glom_images;
$self->rearrange_images();
$self->plot_images;
unlink $self->{POSREG}, $self->{POSCAT};
$self->image_status_pos;
}
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', ''),
);
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} . '1chan', 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,
});
$batbinevt->seriousness(1);
$batbinevt->run();
if( $batbinevt->had_error() == 2 ){
$log->error(2, "batbinevt level 2 error");
}elsif( $batbinevt->had_error() == 1 ){
##################################################################
# 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.
#################################################################
my $errout = $batbinevt->stderr();
foreach ( $errout =~ /^.*$/gm ){
if( ! /^=+$/ && ! /WARNING: no overlapping good time intervals were found/ ){
$log->error(2, "batbinevt level 2 error");
last;
}
}
}
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', ''),
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
pcodemap => 'APPEND_LAST',
);
# 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) {
my $mode = $c->{mode} . $chan;
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
my $name = $filename->get('dpimage', 'bat', 'evpb' . $chan);
$bkgfile = $name if -f $name;
}
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', ''),
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,
});
$sky2xy->seriousness(1);
$sky2xy->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 $imgfits = Util::FITSfile->new($imgfile, 0);
my $xnaxis = $imgfits->keyword('NAXIS1');
my $ynaxis = $imgfits->keyword('NAXIS2');
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;
$xstop = $xnaxis if $xstop > $xnaxis;
$ystop = $ynaxis if $ystop > $ynaxis;
if ($xstart < 1) { $xstart = 1 };
if ($ystart < 1) { $ystart = 1 };
if( $xstop < 1 || $ystop < 1 || $xstart > $xnaxis || $ystart > $ynaxis ){
$log->error(1, "Postage stamp image outside of sky image. " .
"Will not make postage stamp $postfile .");
next;
}
my $copier = Util::HEAdas->new('ftcopy')
->params({
infile => $imgfile . "[$xstart:$xstop,$ystart:$ystop]",
outfile => $postfile,
})
->run;
}
}
}
sub refine_position
{
my ($self) = @_;
my $log = $self->log;
my $filename = $self->filename;
my $procpar = $self->procpar;
my $jobpar = $self->jobpar;
$log->entry('BAT position refinement...');
my $mode = 'evps1chan';
my $ackpath = Util::BATCave::get_tdrss_ack($self);
if (not $ackpath) {
$log->entry('TDRSS ACK unavalable, so no position refinement');
return;
}
my $pointSource = 0;
my $imageTrigger = 0;
my $status = Util::SimpleFITS->readonly($ackpath)
->readkey(POINTSRC => $pointSource)
->readkey(IMAGETRG => $imageTrigger)
->close
->status;
if (not $pointSource or $imageTrigger) {
$log->entry('this is an image trigger or not a point source, so no position refinement');
return;
}
my $imgfile = $filename->get('skyimage', 'bat', $mode, '*');
if (not $imgfile) {
$log->entry('preslew sky image unavalable, so no position refinement');
return;
}
my $pcodefile = $imgfile . '[BAT_PCODE_1]';
my %info = (
ra => $jobpar->read('burst_ra'),
dec => $jobpar->read('burst_dec'),
err_rad => 4 / 60, # deg
incat => 'input.cat',
outcat => $self->{POSCAT},
# $filename->get('srccat', 'bat', '', 0);
region => $self->{POSREG},
imgfile => $imgfile,
pcodefile => $pcodefile,
trigger => $jobpar->read('target'),
);
unlink($info{incat}, $info{outcat});
$self->refine_position_aux(\%info);
}
# BAT::posrefine::refine
sub refine_position_aux
{
my ($self, $info) = @_;
my $log = $self->log;
my $filename = $self->filename;
# Make input catalog
my $incat = $info->{incat};
my $fits = Util::SimpleFITS->create($incat);
my $status = $fits->status();
if ($status) {
$log->error(1, "could not open $incat for writing");
return;
}
# Write catalog columns
$status = $fits
->createtab("BAT_CATALOG")
->writekey("NAXIS2",1)
->insertcol({TTYPE => ["NAME", "Source Name"],
TFORM => "20A"})
->insertcol({TTYPE => ["RA_OBJ", "Source Right Ascension"],
TFORM => "1D", TUNIT => "deg"})
->insertcol({TTYPE => ["DEC_OBJ", "Source Declination"],
TFORM => "1D", TUNIT => "deg"})
->insertcol({TTYPE => ["ERR_RAD", "Source Error Radius"],
TFORM => "1D", TUNIT => "deg"})
->writecol("NAME", {},"TRIG_$info->{trigger}")
->writecol("RA_OBJ", {},$info->{ra})
->writecol("DEC_OBJ",{},$info->{dec})
->writecol("ERR_RAD",{},$info->{err_rad})
->status();
if ($status) {
$log->error(1, "could not create table in $incat");
return;
}
$fits->close();
my $region = $info->{region};
unlink ($region) if ($region ne "NONE" && -f $region);
# Run source detection for this source
my $outcat = $info->{outcat};
my $tool = Util::HEAdas->new('batcelldetect')
->params({
infile => $info->{imgfile},
outfile => $outcat,
snrthresh => 6,
incatalog => $incat,
srcdetect => 'NO',
posfit => 'YES',
distfile => 'CALDB',
niter => 1,
regionfile => $region,
pcodefile => $info->{pcodefile},
clobber => 'yes',
})
->run;
if ($tool->had_error or not -f $outcat) {
$log->error(1, 'batcelldetect failed in position refinement');
return;
}
$fits = Util::SimpleFITS->readonly($outcat)->move("BAT_CATALOG");
$status = $fits->status();
if ($status) {
$log->error(1, "could not open $outcat");
return;
}
my (@ra, @dec, @err_rad, @chi2, @snr);
$status = $fits
->readcol("RA_OBJ", TDOUBLE, [ ], \@ra)
->readcol("DEC_OBJ", TDOUBLE, [ ], \@dec)
->readcol("ERR_RAD", TDOUBLE, [ ], \@err_rad)
->readcol("CHI2_NU", TDOUBLE, [ ], \@chi2)
->readcol("SNR", TDOUBLE, [ ], \@snr)
->status();
$fits->setstatus(0)->close();
if ($status) {
$log->error(1, "could not read $outcat");
return;
}
unless(@ra){
$log->entry("No BAT sources detected.");
}else{
$log->entry("refined position:");
$log->entry(sprintf('RA=%.4f DEC=%.4f', $ra[0], $dec[0]));
$log->entry(sprintf('error=%.4f [arcmin]', $err_rad[0] * 60));
$log->entry(sprintf('chi^2=%.4f', $chi2[0]));
$log->entry(sprintf('snr=%.2f', $snr[0]));
}
###############################################
# Append position refinement to BAT image file
###############################################
Util::FITSfile->new($outcat, 'BAT_CATALOG')
->append_to($info->{imgfile});
unlink $outcat, $incat;
}
sub glomImages
{
my ($self, $info) = @_;
my $log = $self->log;
my $filename = $self->filename;
my $imgtype = $info->{sky} ? 'skyimage' : 'dpimage';
my $img1chan = $filename->get($imgtype, 'bat', $info->{mode} . '1chan');
my $img4chan = $filename->get($imgtype, 'bat', $info->{mode} . '4chan');
if ((not $img1chan or not -f $img1chan)
and (not $img4chan or not -f $img4chan)) {
$log->entry("no $imgtype images");
return;
}
my $glomfile = $filename->get($imgtype, 'bat', $info->{mode}, 0);
$log->entry("packaging $info->{key} ${imgtype}s into $glomfile");
my $appender = Util::HEAdas->new('ftappend')
->params({
outfile => $glomfile,
});
my @error;
if (not $img1chan or not -f $img1chan) {
push(@error, "missing BAT total $imgtype");
rename($img4chan, $glomfile)
or push(@error, "unable to rename $img4chan to $glomfile [$!]");
}
elsif (not $img4chan or not -f $img4chan) {
push(@error, "missing BAT 4 channel $imgtype");
rename($img1chan, $glomfile)
or push(@error, "unable to rename $img1chan to $glomfile [$!]");
}
else {
rename($img4chan, $glomfile)
or push(@error, "unable to rename $img4chan to $glomfile [$!]");
}
foreach my $e (@error) {
$log->error(1, $e);
}
if (not @error) {
my $nhdu = 0;
my $status = Util::SimpleFITS->readonly($img1chan)
->nhdu($nhdu)
->close
->status;
if ($status) {
$log->error(2, "unable to open $img1chan");
}
my $hdupath = $img1chan . '[0]';
if($info->{sky}){
$hdupath .= '[col #EXTNAME = \"BAT_IMAGE_TOT\"]';
}else{
$hdupath .= '[col #EXTNAME = \"BAT_DPI_TOT\"]';
}
$appender->params({
infile => $hdupath,
})->run;
if ($appender->had_error) {
$log->error([ 2, HEATOOL_ERROR ],
"unable to append $hdupath to $glomfile");
}
$status = Util::SimpleFITS->readonly($img1chan)
->move('BAT_CATALOG')
->close
->status;
if (!$status && $info->{sky} && $info->{mode} eq 'evps') {
$hdupath = $img1chan . "[BAT_CATALOG]";
$appender->params({
infile => $hdupath,
})->run;
if ($appender->had_error) {
$log->error([ 2, HEATOOL_ERROR ],
"unable to append $hdupath to $glomfile");
}
}
unlink($img1chan) or $log->error(2, "unable to remove $img1chan [$!]");
}
}
sub glom_images
{
my $self = shift;
my $log = $self->log;
# sky: 1 channel + 4 channel + EBOUNDS + GTI + partial coding map
# det: 1 channel + 4 channel + EBOUNDS + GTI
my @config = (
# DET images
{ key => 'preburst',
mode => 'evpb',
},
{ key => 'preslew',
mode => 'evps',
},
{ key => 'postslew',
mode => 'evas',
},
# SKY images
{ key => 'preburst',
sky => 1,
mode => 'evpb',
},
{ key => 'preslew',
sky => 1,
mode => 'evps',
},
{ key => 'postslew',
sky => 1,
mode => 'evas',
},
);
foreach my $c (@config) {
$self->glomImages($c);
}
}
sub rearrange_images{
my ($self) = @_;
my $log = $self->log;
my $filename = $self->filename;
my $imgtype = 'skyimage';
my @types = qw /evpb evps evas/;
my $stat = Util::HEAdas->new('ftstat');
foreach my $t (@types){
my $img = $filename->get($imgtype, 'bat', $t);
next if (!-e $img);
my $nhdu = 0;
my $status = Util::SimpleFITS->readonly($img)->nhdu($nhdu)->close->status;
my $extnum = undef;
for (my $i=0; $i<$nhdu; $i++){
my $ext = $i;
my $infile = $img."\[$ext\]";
$stat->params({infile => $infile})->run();
my $stdo = $stat->stdout();
my $extname = $self->getExt($stdo, 'extname');
if(defined $extname and $extname =~ /BAT_IMAGE_TOT/){
$extnum = $i;
}
if(defined $extnum){
my $newfile = $img.'_new';
my $fits = Util::FITSfile->new($newfile);
Util::HEAdas->new('ftcopy')
->params({
infile => "$img\[$extnum\]",
outfile => "$newfile",
copyall => 'NO',
})
->run();
Util::HEAdas->new('fthedit')
->params({
infile => $newfile.'[0]',
keyword => 'EXTNAME',
operation => 'add',
value => 'BAT_IMAGE_TOT'
})
->run();
Util::HEAdas->new('fthedit')
->params({
infile => $newfile.'[0]',
keyword => 'HDUNAME',
operation => 'deleteall'
})
->run();
for (my $i=0; $i<$nhdu; $i++){
next if $i == $extnum;
Util::HEAdas->new('ftappend')
->params({
infile => "$img\[$i\]",
outfile => "$newfile",
})
->run();
if($i==0){
Util::HEAdas->new('fthedit')
->params({
infile => $newfile.'[1]',
keyword => 'EXTNAME',
operation => 'add',
value => 'BAT_IMAGE_1'
})
->run();
Util::HEAdas->new('fthedit')
->params({
infile => $newfile.'[1]',
keyword => 'HDUNAME',
operation => 'deleteall'
})
->run();
}
}
if(-e $newfile){
unlink $img;
my $rv = system("mv $newfile $img");
if($rv != 0){
}
}
last;
}
}
}
# exit;
}
sub getExt {
my $self = shift;
my ($str, $pat) = @_;
my @tar = split/\n/, $str;
my $quest = undef;
foreach my $l (@tar){
next if($l =~ /^\s*$/ or $l =~ /^\s*\=/);
$l =~ s/^\s+//;
if($l =~ /$pat\:/){
$quest = (split/\:\s+/, $l)[-1];
last;
}
}
return $quest;
}
sub plot_images {
my ($self) = @_;
my $log = $self->log;
my $filename = $self->filename;
my $jobpar = $self->jobpar;
my $title = uc($jobpar->read('object')) .' SWIFT BAT '.
$jobpar->read('sequence') .' ('.$jobpar->read('obsdate').')';
my $msbce = $filename->get('tdmess', 'bat', 'ce', 0);
my $msbce1 = Util::FITSfile->new($msbce,0);
my $trigdate = undef;
my $trigtime = undef;
$trigdate = $msbce1->keyword('DATETRIG');
$trigtime = $msbce1->keyword('TRIGTIME');
my $trigtimeSrc = "TRIGTIME keyword";
if(!defined $trigtime or $trigtime == 0){
$trigtime = $self->Subs::UvotProduct::getTrigFromDB();
$trigtimeSrc = "TDRSS Catalog";
if(!defined $trigtime){
$trigtime = $self->Subs::UvotProduct::getTrigFromSWIFTCatalog();
$trigtimeSrc = "Lorella Catalog";
if(!defined $trigtime){
$trigtime = $self->Subs::UvotProduct::getTrigFromJDCatalog();
$trigtimeSrc = "JD Catalog";
if(!defined $trigtime){
$trigtime = $self->Subs::UvotProduct::getTrigFromLocalCatalog();
$trigtimeSrc = "Local Catalog";
}
}
}
}
if(!defined $trigdate and defined $trigtime){
my $dateobj = Util::Date->new($trigtime);
$trigdate = $dateobj->date().'T'.$dateobj->time();
}
my $title2 = "Time(s) since BAT trigger time";
if(defined $trigtime){
$title2 .= "(UT $trigdate/MET $trigtime)";
}
my $region = $self->{POSREG};
my @postfiles = $filename->get('postimg', 'bat', '*', 0);
return unless @postfiles;
my %modes = ('evpb' => 'Before Burst', 'evps' => 'Before Slew', 'evas' => 'After Slew');
my $plot_file = $filename->get('skyplot', 'bat', 'ev', 0);;
my $circle;
if( -f $region && open(REG, $region) ){
while(<REG>){
if( /^circle\(([\.\d]+)d,([-\.\d]+)d,([\.\d]+)d\)/ ){
$circle = "ra_dec_to_pixel/ra=$1/dec=$2/color=0/circle=". $3*3600;
}
}
close REG;
}
my $ximage = Util::Xanadu->new("ximage");
my @titlepos = ('/vx=0.2/vy=0.8', '/vx=0.7/vy=0.8', '/vx=0.2/vy=0.4', '/vx=0.7/vy=0.4');
my @titlepos2 = ('/vx=0.12/vy=0.78/csize=0.65', '/vx=0.57/vy=0.78/csize=0.65', '/vx=0.12/vy=0.38/csize=0.65', '/vx=0.57/vy=0.38/csize=0.65');
my $viewports = 'viewports.tmp';
open VIEW, ">$viewports";
print VIEW "0.1 0.45 0.5 0.75\n";
print VIEW "0.55 0.9 0.5 0.75\n";
print VIEW "0.1 0.45 0.1 0.35\n";
print VIEW "0.55 0.9 0.1 0.35\n";
close VIEW;
my @commands;
push @commands, ('cpd /vgif', 'cey 2000');
push @commands, ("viewport/file=$viewports");
push @commands, ('color/setcolor=0/namecolor=white');
push @commands, ('color/setcolor=1/namecolor=black');
push @commands, ('cct/set bb');
push @commands, ('levels/num=239');
my $stat = Util::HEAdas->new('ftstat')
->params({centroid => 'no'});
my $nview = 1;
foreach my $mode (keys %modes){
my $file = ( grep(/${mode}/, @postfiles) )[0];
if($file){
my $fits = Util::FITSfile->new($file,0);
my $Tstart = $fits->keyword('TSTART');
my $Tstop = $fits->keyword('TSTOP');
my $dif1 = sprintf("%.2f", $Tstart-$trigtime);
my $dif2 = sprintf("%.2f", $Tstop-$trigtime);
my $dif3 = sprintf("%.2f", $Tstop - $Tstart);
if($dif1 >= 0.0){
$dif1 = '+'.$dif1;
}
if($dif2 >= 0.0){
$dif2 = '+'.$dif2;
}
# if($dif3 >= 0.0){
# $dif3 = '+'.$dif3;
# }
my $ptitle = 'From T'.$dif1.'s to T'.$dif2.'s = '.$dif3.'s';
$stat->params({infile => $file . '[0]'})
->run();
my $par = $stat->parfile();
push @commands, ("read/fits/mapid=$nview $file");
push @commands, ("disp/lin/mapid=$nview/min=". $par->read('min') ."/max=". $par->read('max'));
push @commands, ("map $nview");
push @commands, ('grid/csize=0.55');
push @commands, ($circle) if $circle;
push @commands, ("label/color=1$titlepos[$nview-1] \"$modes{$mode}\"");
push @commands, ("label/color=1$titlepos2[$nview-1] \"$ptitle\"");
push @commands, ('scale/curvp');
$nview++;
}
}
# push @commands, ("label/vx=0.3/vy=0.9/csize=1.1/color=1 \"BAT Sky Images, ". $jobpar->read('sequence') .'"');
push @commands, ("label/vx=0.17/vy=0.9/csize=1.0/color=1 \"$title\"");
push @commands, ("label/vx=0.125/vy=0.025/csize=0.75/color=1 \"$title2\""); push @commands, 'exit';
$ximage->script(@commands)->run();
rename 'pgplot.gif', $plot_file;
unlink $viewports;
}
sub image_status_pos {
my ($self) = @_;
my $log = $self->log;
my $filename = $self->filename;
my $jobpar = $self->jobpar;
my $tmp = 'img_status.tmp';
my $attitude = $filename->get('attcorr', 'p');
if (-e $attitude) {
$log->entry("Got pat attitude file $attitude");
} else {
$log->error(1, "$attitude pat attitude file does not exist, trying sat file");
$attitude = $filename->get('attitude', 's');
if (-e $attitude) {
$log->entry("Got sat attitude file $attitude");
} else {
$log->error(1, "Unable to find attitude file, exiting");
return;
}
}
my $imgstatpos = Util::HEAdas->new('batimgstatpos')
->params({teldef => 'CALDB',
distfile => 'NONE',
outfile => $tmp,
attitude => $attitude});
foreach my $file ($filename->get('btbimgtr', 'bat', '', '*')){
$imgstatpos->params({infile => $file})
->run();
unless( $imgstatpos->had_error() ){
rename $tmp, $file;
}else{
unlink $tmp;
}
}
}
1;