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.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.8 $ # # ############################################################################## 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 => $filename->fetch_cal('teldef', 'bat'), attfile => $filename->get('attitude', 's'), corrections => 'autocollim,flatfield,ndets,pcode,maskwt', aperture => $filename->fetch_cal('aperture', 'bat'), 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 => $filename->fetch_cal('teldef', 'bat'), attfile => $filename->get('attitude', 's'), # since pcodemap == imageSky(bkgfile=NONE, pcodemap=YES) bkgfile => 'NONE', pcodemap => 'YES', corrections => 'autocollim,flatfield,ndets,pcode,maskwt', aperture => $filename->fetch_cal('aperture', 'bat'), 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;