package Subs::BATSpectra;
##############################################################################
#
# DESCRIPTION: Accumulate BAT spectra and follow on products
#
# Port of Craig Markwardt's BAT::spectra, BAT::lc and setup to SDC.
#
#
# HISTORY: $Log: BATSpectra.pm,v $
# HISTORY: Revision 1.16 2004/12/31 01:36:26 apsop
# HISTORY: Comment out production of spectra which are not yet official.
# HISTORY:
# HISTORY: Revision 1.15 2004/12/06 01:00:12 apsop
# HISTORY: Test for presence of event files, change some error() calls than entry() calls.
# HISTORY:
# HISTORY: Revision 1.14 2004/11/29 22:14:35 apsop
# HISTORY: Corrected efile+escale parameters to batdrmgen. Only process existing
# HISTORY: files.
# HISTORY:
# HISTORY: Revision 1.13 2004/11/16 14:30:25 apsop
# HISTORY: Avoid using undefined value.
# HISTORY:
# HISTORY: Revision 1.12 2004/11/09 22:22:29 apsop
# HISTORY: Corrected name of peak GTI.
# HISTORY:
# HISTORY: Revision 1.11 2004/11/09 22:00:42 apsop
# HISTORY: Collect BAT GTIs in a single file.
# HISTORY:
# HISTORY: Revision 1.10 2004/11/08 18:41:54 apsop
# HISTORY: Implemented Craig's Markwardt's BAT::spectra and BAT::lc modules.
# HISTORY:
#
# VERSION: $Revision: 1.16 $
#
#
##############################################################################
use Subs::Sub;
@ISA = ("Subs::Sub");
use strict;
use Util::HEAdas;
use Util::BATCave;
sub new {
my $proto=shift;
my $self=$proto->SUPER::new();
$self->{DESCRIPTION}= 'Build BAT spectra, response matrices, light curves';
return $self;
}
##################
# METHODS:
##################
sub body {
my $self = shift;
$self->make_spectra;
$self->make_response_matrices;
$self->make_light_curves;
}
sub make_spectra {
my $self = shift;
my $log = $self->log;
my $filename = $self->filename;
my $procpar = $self->procpar;
my $jobpar = $self->jobpar;
$log->entry('Accumulating spectra');
my @events = $filename->getExisting('unfiltered', 'bat', 'evsh{ps,sl,as}');
if (not @events) {
$log->entry('no BAT preslew,slew,postslew events, not creating spectra');
return;
}
# First the full-burst type intervals, use all the events
my %info = (
files => \@events,
qmap => $filename->get('qualcal', 'bat', 'cb'),
);
# $info{ebins} = $filename->fetch_cal('ebounds', 'bat');
#
# $info{outfile} = $filename->get('spectrum', 'bat', 'evto', 0);
# $info{gti} = Util::BATCave::get_gti($self, 'GTI_TOT');
# $self->makeSpectrum(\%info)
# if $info{gti};
#
# $info{outfile} = $filename->get('spectrum', 'bat', 'evt9', 0);
# $info{gti} = Util::BATCave::get_gti($self, 'GTI_T90');
# $self->makeSpectrum(\%info)
# if $info{gti};
#
# $info{outfile} = $filename->get('spectrum', 'bat', 'evt5', 0);
# $info{gti} = Util::BATCave::get_gti($self, 'GTI_T50');
# $self->makeSpectrum(\%info)
# if $info{gti};
#
# $info{outfile} = $filename->get('spectrum', 'bat', 'evpk', 0);
# $info{gti} = Util::BATCave::get_gti($self, 'GTI_PEAK');
# $self->makeSpectrum(\%info)
# if $info{gti};
#
#
###########################################################
# these are Craig's totflu and peakflu spectra...
###########################################################
# $info{ebins} = Util::BATCave::chan4;
#
# $info{outfile} = $filename->get('spectrum', 'bat', 'evtf', 0);
# $info{gti} = Util::BATCave::get_gti($self, 'GTI_TOT');
# $self->makeSpectrum(\%info)
# if $info{gti};
#
# $info{outfile} = $filename->get('spectrum', 'bat', 'evpf', 0);
# $info{gti} = Util::BATCave::get_gti($self, 'GTI_TOT');
# $self->makeSpectrum(\%info)
# if $info{gti};
#
###########################################################
# now the preslew, slew and postslew spectra
###########################################################
$info{ebins} = Util::BATCave::chan4;
$info{outfile} = $filename->get('spectrum', 'bat', 'evps', 0);
$info{gti} = Util::BATCave::get_gti($self, 'GTI_TOT');
$self->makeSpectrum(\%info)
if $info{gti};
$info{outfile} = $filename->get('spectrum', 'bat', 'evsl', 0);
$info{gti} = 'NONE';
$self->makeSpectrum(\%info)
if $info{gti};
$info{outfile} = $filename->get('spectrum', 'bat', 'evas', 0);
$info{gti} = 'NONE';
$self->makeSpectrum(\%info)
if $info{gti};
}
sub makeSpectrum
{
my ($self, $info) = @_;
my $log = $self->log;
my $inlist = join(',', @{ $info->{files} });
my $qmap = $info->{qmap} || 'NONE';
my $gti = $info->{gti} || 'NONE';
unlink($info->{outfile});
##############################################
# set up batbinevt
#############################################
my $binner = Util::HEAdas->new('batbinevt')
->params({
infile => $inlist,
outfile => $info->{outfile},
outtype => 'PHA',
timedel => 0,
timebinalg => 'uniform',
energybins => $info->{ebins},
detmask => $qmap,
ecol => 'ENERGY',
weighted => 'YES',
outunits => 'COUNTS',
gtifile => $info->{gti},
clobber => 'YES',
tstart => 'INDEF',
tstop => 'INDEF',
snrthresh => 6.0,
buffersize => 16384,
chatter => 5,
history => 'YES',
})
->run;
####################
# check for errors
####################
if ($binner->had_error || ! -s $info->{outfile}) {
$log->entry("Did not sucessfully make a spectrum for $inlist");
}
}
sub make_response_matrices
{
my $self = shift;
my $log = $self->log;
my $filename = $self->filename;
my $procpar = $self->procpar;
my $jobpar = $self->jobpar;
$log->entry('Computing response matrices');
my %common = (
calfile => $filename->fetch_cal('rmfparms', 'bat'),
depthfile => $filename->fetch_cal('depthdist', 'bat'),
escale => 'FILE',
efile => $filename->fetch_cal('ebounds', 'bat'),
chatter => 3,
clobber => 'yes',
);
my $batdrmgen = Util::HEAdas->new('batdrmgen')
->params(\%common);
# preslew
my $spectrum = $filename->get('spectrum', 'bat', 'evps');
if (-f $spectrum) {
$batdrmgen->params({
infile => $spectrum,
outfile => $filename->get('response', 'bat', 'evps', 0),
})
->run;
}
# postslew
$spectrum = $filename->get('spectrum', 'bat', 'evas');
if (-f $spectrum) {
$batdrmgen->params({
infile => $spectrum,
outfile => $filename->get('response', 'bat', 'evas', 0),
})
->run;
}
}
sub make_light_curves
{
my $self = shift;
my $log = $self->log;
my $filename = $self->filename;
my $procpar = $self->procpar;
my $jobpar = $self->jobpar;
$log->entry('Accumulating light curves');
my @events = $filename->getExisting('unfiltered', 'bat', 'evsh{ps,sl,as}');
if (not @events) {
$log->entry('no BAT preslew,slew,postslew events, not creating light curves');
return;
}
my @config = (
{ desc => 'weighted 4ms light curve',
outfile => $filename->get('lightcurve', 'bat', 'evbu', 0),
ebins => Util::BATCave::chan1,
gti => 'NONE',
tbinsize => 0.004,
weighted => 'YES',
},
{ desc => 'weighted 64ms 1-channel light curve',
outfile => $filename->get('lightcurve', 'bat', 'evms', 0),
ebins => Util::BATCave::chan1,
gti => 'NONE',
tbinsize => 0.064,
weighted => 'YES',
},
{ desc => 'unweighted 4ms light curve',
outfile => $filename->get('rawlc', 'bat', 'evbu', 0),
ebins => Util::BATCave::chan1,
gti => 'NONE',
tbinsize => 0.004,
weighted => 'NO',
},
);
my $bbfile = Util::BATCave::get_gti($self, 'GTI_BAYES');
push(@config,
{ desc => 'weighted 4-channel light curve, Bayesian blocks',
outfile => $filename->get('lightcurve', 'bat', 'evbb4c', 0),
ebins => Util::BATCave::chan4,
gti => $bbfile,
tbinsize => 'g',
weighted => 'YES',
},
{ desc => 'weighted 1-channel light curve, Bayesian blocks',
outfile => $filename->get('lightcurve', 'bat', 'evbb1c', 0),
ebins => Util::BATCave::chan1,
gti => $bbfile,
tbinsize => 'g',
weighted => 'YES',
},
) if ($bbfile and -f $bbfile);
my $qmap = $filename->get('qualcal', 'bat', 'cb') || 'NONE';
foreach my $c (@config) {
$c->{files} = \@events;
$c->{qmap} = $qmap;
$self->lcMake($c);
}
}
# BAT::lc->make
sub lcMake
{
my ($self, $info) = @_;
my $log = $self->log;
my $filename = $self->filename;
my $procpar = $self->procpar;
my $jobpar = $self->jobpar;
my %par = (
infile => join(',', @{ $info->{files} }),
outfile => $info->{outfile},
outtype => 'LC',
energybins => $info->{ebins},
detmask => $info->{qmap},
weighted => $info->{weighted} || 'NO',
outunits => 'RATE',
ecol => 'ENERGY',
gtifile => $info->{gti} || 'NONE',
clobber => 'yes',
);
if ($info->{tbinsize} =~ /^g/i) {
$par{timedel} = 0;
$par{timebinalg} = 'gti';
}
else {
$par{timedel} = $info->{tbinsize};
$par{timebinalg} = 'uniform';
}
my $tool = Util::HEAdas->new('batbinevt')
->params(\%par)
->run;
}